{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "qaoa_hands_on_pubV2.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/JackHidary/quantumcomputingbook/blob/master/Notes/qaoa_hands_on_pubV2.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gAxVvd1LF3l4",
        "colab_type": "text"
      },
      "source": [
        "# Quantum Approximate Optimization Algorithm\n",
        "\n",
        "This notebook provides an introduction to the Quantum Approximate Optimization Algorithm (QAOA) using the Cirq. The presentation mostly follows [Farhi et al.](https://arxiv.org/abs/1411.4028). We will show how to construct the QAOA circuit and use it to solve some simple problems.\n",
        "\n",
        "## Instructions\n",
        "\n",
        "Go to File --> Save Copy in Drive to get your own copy to play with."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "RlJBDvNgC00H",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "# install cirq\n",
        "!pip install cirq==0.5 --quiet"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GPjUqrK8DJTq",
        "colab_type": "text"
      },
      "source": [
        "To verify that Cirq is installed in your environment, try to `import cirq` and print out a diagram of the Bristlecone device."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "FTrmLyq4C2gf",
        "colab_type": "code",
        "outputId": "48bc4f5f-6575-4aef-dc20-ab167f04ef7a",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 575
        }
      },
      "source": [
        "import cirq\n",
        "import numpy as np\n",
        "import sympy\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "print(cirq.google.Bristlecone)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "                                             (0, 5)────(0, 6)\n",
            "                                             │         │\n",
            "                                             │         │\n",
            "                                    (1, 4)───(1, 5)────(1, 6)────(1, 7)\n",
            "                                    │        │         │         │\n",
            "                                    │        │         │         │\n",
            "                           (2, 3)───(2, 4)───(2, 5)────(2, 6)────(2, 7)───(2, 8)\n",
            "                           │        │        │         │         │        │\n",
            "                           │        │        │         │         │        │\n",
            "                  (3, 2)───(3, 3)───(3, 4)───(3, 5)────(3, 6)────(3, 7)───(3, 8)───(3, 9)\n",
            "                  │        │        │        │         │         │        │        │\n",
            "                  │        │        │        │         │         │        │        │\n",
            "         (4, 1)───(4, 2)───(4, 3)───(4, 4)───(4, 5)────(4, 6)────(4, 7)───(4, 8)───(4, 9)───(4, 10)\n",
            "         │        │        │        │        │         │         │        │        │        │\n",
            "         │        │        │        │        │         │         │        │        │        │\n",
            "(5, 0)───(5, 1)───(5, 2)───(5, 3)───(5, 4)───(5, 5)────(5, 6)────(5, 7)───(5, 8)───(5, 9)───(5, 10)───(5, 11)\n",
            "         │        │        │        │        │         │         │        │        │        │\n",
            "         │        │        │        │        │         │         │        │        │        │\n",
            "         (6, 1)───(6, 2)───(6, 3)───(6, 4)───(6, 5)────(6, 6)────(6, 7)───(6, 8)───(6, 9)───(6, 10)\n",
            "                  │        │        │        │         │         │        │        │\n",
            "                  │        │        │        │         │         │        │        │\n",
            "                  (7, 2)───(7, 3)───(7, 4)───(7, 5)────(7, 6)────(7, 7)───(7, 8)───(7, 9)\n",
            "                           │        │        │         │         │        │\n",
            "                           │        │        │         │         │        │\n",
            "                           (8, 3)───(8, 4)───(8, 5)────(8, 6)────(8, 7)───(8, 8)\n",
            "                                    │        │         │         │\n",
            "                                    │        │         │         │\n",
            "                                    (9, 4)───(9, 5)────(9, 6)────(9, 7)\n",
            "                                             │         │\n",
            "                                             │         │\n",
            "                                             (10, 5)───(10, 6)\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "3wxPYfCo6Ygw",
        "colab_type": "text"
      },
      "source": [
        "### Description of QAOA\n",
        "\n",
        "Let's start out with a description of the QAOA algorithm. We'll discuss the structure of the problem it tries to solve and the quantum circuit we need to build to implement it."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "8A7a3jcql1l5",
        "colab_type": "text"
      },
      "source": [
        "Suppose you have a function $C(z)$ that depends on a collection of variables $z = z_1,z_2,\\ldots, z_n$, where each $z_j$ can be equal to either $+1$ or $-1$ (the important thing here is that each $z_j$ has two possible values, and by convention we choose those values to be $\\pm 1$). The QAOA is a general-purpose algorithm whose goal is to produce an assignment of the $z_j$ that gives a relatively low value of $C(z)$. It's not guaranteed to give the lowest possible value of $C(z)$---hence the name \"Approximate\"---except in a particular limit which we will discuss.\n",
        "\n",
        "The QAOA algorithm acts on $n$ qubits. As you might guess, each qubit represents one of the variables in our function, and the $2^n$ states of the computational basis correspond to the $2^n$ possible assignments of the $z$ variables. To be more specific, let's agree that the value of $z_j$ corresponds to the measurement outcome of the Pauli-$Z$ operator on the $j$th qubit. There is a potential confusion here because the state $| 0 \\rangle$ corresponds to $z = +1$, while the state $| 1\\rangle$ corresponds to $z=-1$. This is unfortunate, but is something that we'll just have to deal with."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "q4F7ccpiGOn4",
        "colab_type": "text"
      },
      "source": [
        "The QAOA algorithm is fairly simple to explain, though the reasons behind why it works are not obvious at first glance. As usual, we begin with all of our qubits initialized in the $|0\\rangle$ state. The first step is to act with $H^{\\otimes n}$, the Hadamard operator on each qubit. This prepares an equal superposition of all bitstrings, i.e., an equal superposition of all possible $z$ assignments:\n",
        "$$\n",
        "H^{\\otimes n} |0^n\\rangle =\\frac{1}{2^{n/2}} \\sum_{x \\in \\{0,1\\}^n} |x\\rangle.\n",
        "$$\n",
        "This should be thought of as the \"real\" initial state of the algorithm. The point of the remaining steps is to affect the amplitudes so that those with small $C(z)$ values grow while those with large $C(z)$ values shrink. Then at the end when we measure the qubits we'll be more likely to find a bitstring with a small value of $C(z)$."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AW1bxjmdHw44",
        "colab_type": "text"
      },
      "source": [
        "The meat of the algorithm relies on the following unitary operator:\n",
        "$$\n",
        "U(\\gamma,C) = e^{i\\pi \\gamma C(Z)/2}.\n",
        "$$\n",
        "This operator deserves some explanation. First, $\\gamma$ is a paramter which we will later treat as a variational parameter, adjusting its value to produce the best possible result. $C$ here is the function we are trying to minimize, and the notation $C(Z)$ is supposed to tell you to plug in the Pauli-$Z$ operator for each qubit in place of the argument $z$. For example, if \n",
        "$$\n",
        "C(z) = 3z_1 z_2 - z_2z_3 + z_4\n",
        "$$\n",
        "then\n",
        "$$\n",
        "C(Z) = 3Z_1 Z_2 - Z_2Z_3 + Z_4.\n",
        "$$\n",
        "It looks like I didn't do much, but the point here is that $C(z)$ is a number while $C(Z)$ is a matrix. That matrix is diagonal in the computational basis, and those diagonal entries represent all the possible values of $C(z)$.\n",
        "\n",
        "After acting with $H^{\\otimes n}$, we act with $U(C,\\gamma)$. The result is still a sum over all possible bit-strings, but now the coefficients are complex phases which depend on $C$. At this point there is still an equal probability to measure any particular string, though, because Born's rule only depends on the square of the amplitude. So the algorithm is not done yet. Below we will have to figure out how to implement $U(\\gamma, C)$ in Cirq so that we can perform this step in the algorithm."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ocxry8MKJFeD",
        "colab_type": "text"
      },
      "source": [
        "The next step of the algorithm is to act with the unitary operator\n",
        "$$\n",
        "U(\\beta,B) = e^{i\\pi\\beta B/2},~~~ B = \\sum_{j=1}^n X_j,\n",
        "$$\n",
        "where $\\beta$ is another variational parameter. Since the Pauli-$X$ operators on each qubit commute with each other, we can alternatively write this as\n",
        "$$\n",
        "U(\\beta, B) = \\prod_{j=1}^n e^{i\\pi\\beta X_j/2}.\n",
        "$$\n",
        "So this is just a rotation of each qubit around the $X$-axis on the Bloch sphere by an amount determined by $\\beta$. This operation is _not_ diagonal in the computational basis, and the resulting state will not be an equal superposition over all bitstrings. So after this step there will be constructive and destructive interference, which hopefully leads to enhancement of states corresponding to small values of $C$. This $U(\\beta, B)$ is sometimes called a \"mixing\" operation. Note that, up to an inconsequential global phase, we can also write\n",
        "$$\n",
        "U(\\beta, B) = \\prod_{j=1}^n X_j^{\\beta}.\n",
        "$$"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "RFgdxbSlLh0y",
        "colab_type": "text"
      },
      "source": [
        "The total circuit consists of repeating the previous two steps a total of $p\\geq 1$ times, where the choice of $p$ is up to you. The parameters $\\gamma$ and $\\beta$ can be chosen independently at each step. So at the conclusion of the circuit, the state of the qubits is\n",
        "$$\n",
        "|\\gamma,\\beta\\rangle = U(\\beta_p,B)U(\\gamma_p,C)\\cdots U(\\beta_1,B)U(\\gamma_1,C)H^{\\otimes n}|0^n\\rangle.\n",
        "$$\n",
        "If we choose $\\gamma$ and $\\beta$ so that the expectation value\n",
        "$$\n",
        "F(\\gamma,\\beta) = \\langle \\gamma,\\beta|C(Z)|\\gamma,\\beta\\rangle\n",
        "$$\n",
        "is minimized, then measuring the state $|\\gamma,\\beta\\rangle$ in the computational basis gives us a good candidate bitstring for the minimum of $C(z)$. That's the whole thing!"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MK02MlydMUqL",
        "colab_type": "text"
      },
      "source": [
        "In summary we have to perform the following tasks in order to implement the QAOA:\n",
        "\n",
        "\n",
        "1.   Figure out out to perform the $U(\\gamma, C)$ operation in Cirq for our choice of $C$.\n",
        "2.   Create a quantum circuit alternating $U(\\gamma, C)$ and $U(\\beta, B)$ operations as many times as desired. \n",
        "3.   Find the optimal value of the variational parameters in our circuit.\n",
        "4.   Measure the output of our circuit."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "X15yPl_KQ20Z",
        "colab_type": "text"
      },
      "source": [
        "### Toy problem: ground state of Ising model\n",
        "\n",
        "The Ising Model defines the energy function\n",
        "$$\n",
        "E = -\\sum_{\\langle i,j \\rangle} Z_i Z_j - \\sum_i h_i Z_i,\n",
        "$$\n",
        "where the notation $\\langle i,j\\rangle$ means a sum over all nearest-neighbor pairs. The picture here is that the qubits live on the vertices of a graph, and the edges of the graph define which qubits are neighbors. We'll just take out graph to be a rectangular lattice with some number of rows and some number of columns. The numbers $h_i$ have the physical interpretation of an external magnetic field.\n",
        "\n",
        "We are interested in finding a low-lying state of the Ising Model, by which I mean a state that has a relatively low amount of energy. This is a difficult problem in general. The pairwise interaction terms would tell you that neighboring qubits should be in the same state to lower the energy, while the magnetic field terms tell you that a given qubit wants to point \"in the same direction\" as its local field, and the strength of that preference depends on the magnitude of the field. These two different kinds of pressure are not always in agreement!\n",
        "\n",
        "This type of problem is a perfect candidate for the QAOA, where we use the energy $E$ as our cost function $C$."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "VFwmWPf7D057",
        "colab_type": "text"
      },
      "source": [
        "### ZZ Gate\n",
        "\n",
        "The first thing we need to do is create the operation $U(\\gamma, C)$, where $C$ is equal to the Ising Model energy function. The first thing to note is that, since all of the terms in the energy commute, we can decompose this operation as\n",
        "$$\n",
        "U(\\gamma, C) = \\prod_{\\langle i,j\\rangle}e^{-i\\pi\\gamma Z_iZ_j/2} \\prod_i e^{-i\\pi \\gamma h_i Z_i/2}.\n",
        "$$\n",
        "This requires that we have the two-qubit gate $\\exp(-i\\pi\\gamma ZZ/2)$ at our disposal. In matrix form, this is\n",
        "$$\n",
        "\\begin{align}\n",
        "\\exp(-i \\pi\\gamma Z\\otimes Z/2) = \\begin{bmatrix}\n",
        "e^{-i\\pi \\gamma/2}  & 0  &0&0\\\\\n",
        "0 & e^{i\\pi \\gamma/2}   &0&0\\\\\n",
        "0&0& e^{i\\pi \\gamma/2}  &0 \\\\\n",
        "0&0 & 0  & e^{-i\\pi \\gamma/2}\n",
        "\\end{bmatrix}\n",
        "\\end{align}\n",
        "$$\n",
        "As of version 0.5.0, Cirq has a built-in gate `cirq.ZZ` which is equivalent to this once you account for a global phase:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3CZsDwTlDy_h",
        "colab_type": "code",
        "outputId": "6626bcf2-f68f-4b1c-c4e0-e716dc8e4f76",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 143
        }
      },
      "source": [
        "a = cirq.NamedQubit(\"a\")\n",
        "b = cirq.NamedQubit(\"b\")\n",
        "gamma = 0.3 # Put your own value here.\n",
        "circuit = cirq.Circuit.from_ops(cirq.ZZ(a,b)**gamma)\n",
        "print(circuit)\n",
        "cirq.unitary(circuit).round(2)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "a: ───ZZ───────\n",
            "      │\n",
            "b: ───ZZ^0.3───\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "array([[1.  +0.j  , 0.  +0.j  , 0.  +0.j  , 0.  +0.j  ],\n",
              "       [0.  +0.j  , 0.59+0.81j, 0.  +0.j  , 0.  +0.j  ],\n",
              "       [0.  +0.j  , 0.  +0.j  , 0.59+0.81j, 0.  +0.j  ],\n",
              "       [0.  +0.j  , 0.  +0.j  , 0.  +0.j  , 1.  +0.j  ]])"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 3
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "WvFVcT2de0V6",
        "colab_type": "text"
      },
      "source": [
        "We should also check that the matrix is what we expect:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "CazWJFM5eGgB",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "test_matrix = np.array([[np.exp(-1j*np.pi*gamma/2),0, 0, 0],\n",
        "                        [0, np.exp(1j*np.pi*gamma/2), 0, 0],\n",
        "                        [0, 0, np.exp(1j*np.pi*gamma/2), 0],\n",
        "                        [0, 0, 0,np.exp(-1j*np.pi*gamma/2)]])\n",
        "cirq.testing.assert_allclose_up_to_global_phase(test_matrix, cirq.unitary(circuit), atol=1e-5)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "cy1ns-RLMzp3",
        "colab_type": "text"
      },
      "source": [
        "### Z Gate\n",
        "\n",
        "The magnetic field terms can be handled in a similar way. The single-qubit unitary\n",
        "$$\n",
        "\\exp(-i\\pi \\gamma hZ/2) = \\begin{bmatrix}\n",
        "e^{-i\\pi \\gamma h/2} & 0 \\\\\n",
        "0 & e^{i\\pi \\gamma h/2}\n",
        "\\end{bmatrix}\n",
        "$$\n",
        "is equivalent up to global phase to `cirq.Z**(h*gamma)`:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "r4BEICATObhU",
        "colab_type": "code",
        "outputId": "7dbe972c-9623-4ba0-af8c-518823301887",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 71
        }
      },
      "source": [
        "a = cirq.NamedQubit(\"a\")\n",
        "gamma = 0.3 # Put your own value here.\n",
        "h = 1.3 # Put your own value here.\n",
        "circuit = cirq.Circuit.from_ops(cirq.Z(a)**(gamma*h))\n",
        "print(circuit)\n",
        "print(cirq.unitary(circuit).round(2))\n",
        "\n",
        "test_matrix = np.array([[np.exp(-1j*np.pi*gamma*h/2),0],\n",
        "                        [0, np.exp(1j*np.pi*gamma*h/2)]])\n",
        "cirq.testing.assert_allclose_up_to_global_phase(test_matrix, cirq.unitary(circuit), atol=1e-5)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "a: ───Z^0.39───\n",
            "[[1.  +0.j   0.  +0.j  ]\n",
            " [0.  +0.j   0.34+0.94j]]\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "esgOuHF_FT_2",
        "colab_type": "text"
      },
      "source": [
        "### Exercise: More general two-qubit gate\n",
        "\n",
        "The Ising Model is particularly simple because the nearest-neighbor interaction $Z_i Z_j$ is already given in terms of a product of Pauli matrices. But suppose instead I told you that the cost function was a sum of terms that looked like\n",
        "$$\n",
        "C(z_i,z_j) = \\begin{cases}\n",
        "c_{00} \\text{ if } z_i =1,~z_j=1,\\\\\n",
        "c_{01} \\text{ if } z_i =1,~z_j=-1,\\\\\n",
        "c_{10} \\text{ if } z_i =-1,~z_j=1,\\\\\n",
        "c_{11} \\text{ if } z_i =-1,~z_j=-1.\n",
        "\\end{cases}\n",
        "$$\n",
        "For some numbers $c_{ab}$. How would you make the analogous two-qubit gate for this case?\n",
        "\n",
        "You can either make a custom gate from scratch, or build a solution from the standard elementary gates."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uaDb6B_jPgrb",
        "colab_type": "text"
      },
      "source": [
        "### Create Circuit\n",
        "\n",
        "Using the `cirq.ZZ` gate we can now create the QAOA circuit. We're going to focus on the Ising Model an a rectangular lattice with an arbitrary number of rows and columns. Here are some things to think about:\n",
        "\n",
        "1.   `cirq.GridQubit`s are natural because our qubits actually do live on a grid. Cirq does not care what kind of qubit you make, though.\n",
        "2.   It's a good idea to define separate functions to place the C and B layers for the circuit. Really these should be generators that yield the required gates.\n",
        "3.   You might consider wrapping everything inside a class. We won't do that here, but if you want to play around with different numbers of rows/columns or different numbers of B/C layers it can be convenient.\n",
        "\n",
        "First we'll define the basic paramters of our model and the generators for the different layers:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "NPeFsF6yFvGr",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "n_cols = 3\n",
        "n_rows = 3\n",
        "h = 0.5*np.ones((n_rows,n_cols))\n",
        "\n",
        "# Arranging the qubits in a list-of-lists like this makes them easy to refer to later.\n",
        "qubits = [[cirq.GridQubit(i,j) for j in range(n_cols)] for i in range(n_rows)]\n",
        "\n",
        "\n",
        "def beta_layer(beta):\n",
        "  \"\"\"Generator for U(beta, B) layer (mixing layer) of QAOA\"\"\"\n",
        "  for row in qubits:\n",
        "    for qubit in row:\n",
        "      yield cirq.X(qubit)**beta\n",
        "    \n",
        "def gamma_layer(gamma, h):\n",
        "  \"\"\"Generator for U(gamma, C) layer of QAOA\n",
        "  \n",
        "  Args:\n",
        "    gamma: Float variational parameter for the circuit\n",
        "    h: Array of floats of external magnetic field values\n",
        "  \"\"\"\n",
        "  for i in range(n_rows):\n",
        "    for j in range(n_cols):\n",
        "      if i < n_rows-1:\n",
        "        yield cirq.ZZ(qubits[i][j], qubits[i+1][j])**gamma\n",
        "      if j < n_cols-1:\n",
        "        yield cirq.ZZ(qubits[i][j], qubits[i][j+1])**gamma\n",
        "      yield cirq.Z(qubits[i][j])**(gamma*h[i,j])"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "M-3Tg-_dlOn0",
        "colab_type": "text"
      },
      "source": [
        "Let's test these functions by constructing the circuit. Try making a circuit with different numbers of layers. How would you automatically make a circuit with a specified number of layers? Make sure the parameters of these layers are distinct `sympy.Symbol`s for later optimization. Print the circuit to see that it's doing what you want it to do."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "4ZmjE6pmtUoG",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "qaoa = cirq.Circuit()\n",
        "qaoa.append(cirq.H.on_each(*[q for row in qubits for q in row]))\n",
        "# YOUR CODE HERE\n",
        "print(qaoa)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "VEAt5QZvtPu_",
        "colab_type": "text"
      },
      "source": [
        "#### Solution"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7zWHPT1ktlUk",
        "colab_type": "text"
      },
      "source": [
        "We'll just illustrate the solution for a single $C$ layer and a single $B$ layer."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "lHjIRxL13nXP",
        "colab_type": "code",
        "outputId": "88ec3ae0-13da-4d2c-8aa0-95e8b2909921",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 379
        }
      },
      "source": [
        "qaoa = cirq.Circuit()\n",
        "gamma = sympy.Symbol('g')\n",
        "beta = sympy.Symbol('b')\n",
        "qaoa.append(cirq.H.on_each(*[q for row in qubits for q in row]))\n",
        "qaoa.append(gamma_layer(gamma,h))\n",
        "qaoa.append(beta_layer(beta))\n",
        "print(qaoa)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "                                                ┌──────────────────┐               ┌──────────────────┐\n",
            "(0, 0): ───H───ZZ─────ZZ─────Z^(0.5*g)───X^b────────────────────────────────────────────────────────────────────────────────────────────\n",
            "               │      │\n",
            "(0, 1): ───H───┼──────ZZ^g───ZZ──────────ZZ──────Z^(0.5*g)─────────────X^b──────────────────────────────────────────────────────────────\n",
            "               │             │           │\n",
            "(0, 2): ───H───┼─────────────┼───────────ZZ^g────ZZ────────────────────Z^(0.5*g)────X^b─────────────────────────────────────────────────\n",
            "               │             │                   │\n",
            "(1, 0): ───H───ZZ^g───ZZ─────┼───────────ZZ──────┼────────Z^(0.5*g)────X^b──────────────────────────────────────────────────────────────\n",
            "                      │      │           │       │\n",
            "(1, 1): ───H──────────┼──────ZZ^g────────ZZ^g────┼────────ZZ───────────ZZ───────────Z^(0.5*g)─────────────X^b───────────────────────────\n",
            "                      │                          │        │            │\n",
            "(1, 2): ───H──────────┼──────────────────────────ZZ^g─────┼────────────ZZ^g─────────ZZ────────────────────Z^(0.5*g)───X^b───────────────\n",
            "                      │                                   │                         │\n",
            "(2, 0): ───H──────────ZZ^g────────────────────────────────┼────────────ZZ───────────┼────────Z^(0.5*g)────X^b───────────────────────────\n",
            "                                                          │            │            │\n",
            "(2, 1): ───H──────────────────────────────────────────────ZZ^g─────────ZZ^g─────────┼─────────────────────ZZ──────────Z^(0.5*g)───X^b───\n",
            "                                                                                    │                     │\n",
            "(2, 2): ───H────────────────────────────────────────────────────────────────────────ZZ^g──────────────────ZZ^g────────Z^(0.5*g)───X^b───\n",
            "                                                └──────────────────┘               └──────────────────┘\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "3HtlMxa6QpVo",
        "colab_type": "text"
      },
      "source": [
        "### Define Expectation Value\n",
        "\n",
        "To train the QAOA circuit---that is, find the optimal values of the paramters---we're going to need to be able to compute the expectation value of the Ising Model energy. We'll do this within Cirq by defining an energy function. We'll divide the total energy by the number of qubits to keep the numbers under control, basically because we expect the energy to scale with the size of the system.\n",
        "\n",
        "If we were using real hardware, the only way to compute the expectation value of the energy would be to estimate it by sampling. Using the simulator we can alternatively compute the wavefunction and then get calculate the expectation value from that. Not only does this save us from having to worry about statistical error, it also tends to be faster that simulating the sampling process."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab_type": "code",
        "id": "-9etj1AeK6dG",
        "colab": {}
      },
      "source": [
        "def energy_from_wavefunction(wf, h):\n",
        "  \"\"\"Computes the energy-per-site of the Ising Model directly from the\n",
        "  a given wavefunction. \n",
        "  \n",
        "  Args:\n",
        "    wf: Array of size 2**(n_rows*n_cols) specifying the wavefunction.\n",
        "    h: Array of shape (n_rows, n_cols) giving the magnetic field values.\n",
        "  \n",
        "  Returns:\n",
        "    energy: Float equal to the expectation value of the energy per site \n",
        "  \"\"\"\n",
        "  \n",
        "  n_sites = n_rows*n_cols\n",
        "\n",
        "  # Z is an array of shape (n_sites, 2**n_sites). Each row consists of the \n",
        "  # 2**n_sites non-zero entries in the operator that is the Pauli-Z matrix on\n",
        "  # one of the qubits times the identites on the other qubits. The\n",
        "  # (i*n_cols + j)th row corresponds to qubit (i,j).\n",
        "  Z = np.array([(-1)**(np.arange(2**n_sites) >> i) for i in range(n_sites-1,-1,-1)])\n",
        "\n",
        "  # Create the operator corresponding to the interaction energy summed over all\n",
        "  # nearest-neighbor pairs of qubits\n",
        "  ZZ_filter = np.zeros_like(wf, dtype=float)\n",
        "  for i in range(n_rows):\n",
        "    for j in range(n_cols):\n",
        "      if i < n_rows-1:\n",
        "        ZZ_filter += Z[i*n_cols + j]*Z[(i+1)*n_cols + j]\n",
        "      if j < n_cols-1:\n",
        "        ZZ_filter += Z[i*n_cols + j]*Z[i*n_cols + (j+1)]\n",
        "\n",
        "  energy_operator = -ZZ_filter - h.reshape(n_sites).dot(Z)\n",
        "        \n",
        "        \n",
        "  # Expectation value of the energy divided by the number of sites\n",
        "  return np.sum(np.abs(wf)**2 * energy_operator) / n_sites"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fjFPEQuyvxjR",
        "colab_type": "text"
      },
      "source": [
        "We'll also need a helper function that computes the expected value of the energy given some parameters of the QAOA."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "XOYLY_u5K7z0",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def energy_from_params(gamma, beta, qaoa, h):\n",
        "  sim = cirq.Simulator()\n",
        "  params = cirq.ParamResolver({'g':gamma, 'b':beta})\n",
        "  wf = sim.simulate(qaoa, param_resolver = params).final_state\n",
        "  return energy_from_wavefunction(wf, h)  "
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "r-CjbPwkRI_I",
        "colab_type": "text"
      },
      "source": [
        "### Training\n",
        "\n",
        "Now we need to figure out the best values of $\\gamma$ and $\\beta$ by minimizing the expectation value of the energy. We'll start by doing a brute-force search of the parameter space for illustrative purposes:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "hM2Zd_kTI578",
        "colab_type": "code",
        "outputId": "701a307e-8aa4-450e-c029-49840a3a7028",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 53
        }
      },
      "source": [
        "%%time\n",
        "grid_size = 50\n",
        "gamma_max = 2\n",
        "beta_max = 2\n",
        "\n",
        "energies = np.zeros((grid_size,grid_size))\n",
        "for i in range(grid_size):\n",
        "  for j in range(grid_size):\n",
        "    energies[i,j] = energy_from_params(i*gamma_max/grid_size, j*beta_max/grid_size, qaoa, h)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "CPU times: user 28.6 s, sys: 15.9 ms, total: 28.6 s\n",
            "Wall time: 28.6 s\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AFP2Ofi0KTfq",
        "colab_type": "code",
        "outputId": "3e579d83-e394-47c8-f088-9f9eba93dc3f",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        }
      },
      "source": [
        "plt.ylabel('gamma')\n",
        "plt.xlabel('beta')\n",
        "plt.title('Energy as a Function of Parameters')\n",
        "plt.imshow(energies, extent = (0,beta_max,gamma_max,0));"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAARsAAAEWCAYAAABWszP/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAIABJREFUeJztnXm8LUV177+/vc94B/DCFRVE0IiJ\naBSR4JiILwbRoJhonjhjNESN0ZeYZzSDGof3jEk0AxpFJahxjEZFH6gkiiYqChgcwKgIGkCQ4TLc\n+Z6z93p/dJ+b3avqnO6779m9zzms7+dzP3d3d1VXdVWf2lW/vVYtmRlBEASjpjPuCgRBcMcgBpsg\nCFohBpsgCFohBpsgCFohBpsgCFohBpsgCFohBptgrEh6hqTPjaHcR0j6gaRtkp7Udvl3RNbcYCPp\nR5J2li/Rwr8zxl2vlYiksyXtcW311BGWd6QkkzSxcM7M3m9mJ46qzCV4LXCGmW0ws0/4i+49+mnZ\nVhvGUM99oqzn68ddjxxrbrApeUL5Ei38e/FyFzD4B7PKeZNrqw+Pu0ItcQRwWU2aJ5jZBuBY4Djg\nT/a1kNX2noyyvmt1sMki6TRJ/y7pLyXdIukqSY8buH6gpHdLuk7StZJeL6k7kPfLkt4i6WbgNZK6\nkv5K0k3lvV688M0t6TckXeLK/31Jn1ykbs+V9F1JWyVdKem3B65tlvRpSbdK2iLp3yRl+07S30i6\nWtLtki6R9ItDtpVJuvfA8d5vTEknSLpG0ssk3VC213MH0s6W7fJjSbeVbT4LfKlMcms5Y3jYQp8M\n5H24pIvKfBdJevjAtQskva7sh62SPidp8xLP8FuSrijb7BxJh5bnfwjcC/hUWY/ppdrCzK4FzgPu\nX+Zfqq8W2uYPJV0P/IOkTWX/3Vi+d5+WdHf3XK+X9JWyPp+SdLCk95f9eJGkIwfS/5yk88vn+p6k\n/1mePx14BvDyhfuU5w+V9LGy/KskvWTgXq+R9FFJ/yjpduA0ScdLurgs+6eS3rxU+zTGzNbUP+BH\nwGMWuXYaMAf8FtAFXgj8BFB5/ePAO4D1wCHA14HfHsg7D/wuMAHMAi8ALgfuDmwC/gWw8vo0sAW4\n70D5/wE8eZG6/SrwM4CARwE7gGPLa/8XeDswWf77xYU6Z+7zTODgsg4vA64HZhZJezbw+kWuGXDv\nXFrghLItXlvW5/FlfTeV198KXAAcVrbzw8v2OHKhfVyf/Hv5+SDgFuBZZf2fVh4fXF6/APghcJ+y\n/S8A3rhI/f8HcBPFrGQa+DvgS03eE38dOJxiFvS6Bn210DZ/XpY7W/bHk4F1wEbgn4BPDJR1AXBF\nec8DKd6p7wOPKdvhvcA/lGnXA1cDzy2vPah8zqNzfUoxobgEeBUwRTHIXgk8trz+Goq/iSeVaWeB\nrwLPKq9vAB66LH+b4x4clvtf+ZJsA24d+PdbAy/2FQNp15Uv/12BuwC7gdmB608DvjCQ979cWZ+n\nHIzK48cw8McE/D3whvLz/Sj+cKYbPscngJeWn18LfJKBP/59aI9bgAcucu1sYNdAO900cK1usNlJ\nddC4AXho+cLuzJVJ/WDzLODrLs9XgdPKzxcAfzJw7UXAZxZ5tndTLBEXjjeUf1RHDrwndYPNwnv0\nY+Btg+/GEn11ArCHRQb4Ms0xwC0DxxcAfzxw/FfAeQPHTwAuLT8/Ffg3d793AK/2/VQePyTz3r6S\n/x68XsPAIFye+xLwZ8Dm5fzbXKvLqCeZ2Z0G/r1z4Nr1Cx/MbEf5cQPFGn4SuK5crtxK0YmHDOS9\n2pVzqDvnr78HeLokUfwhfcTMducqLOlxki4sp8a3UswWFpYIf0Hxzfe5ctr+isUeXNIflFP828r7\nHDhwnxx/OdBOS6Xz3Gxm8wPHOyjacTMwQzED2VcOpfjDHuTHFDOkBa4f+LxQZu29zGwbcLO7Vx0L\n79ERZvYiM9sJtX0FcKOZ7Vo4kLRO0jvKZeXtFH/Md1K5RC/56cDnnZnjhec8AnjIwjtalv8Mii/M\nHEcAh7r0f0Tx5bqAf2+fRzF7/M9yCXfyYg20L6wq8WrEXE0xs9ns/ogG8S7y11EsoRY4vJLY7EJJ\neyiWPU8v/yWUmsHHgGcDnzSzOUmfoJimY2ZbKZZEL5N0f+Dzki4ys3919/lF4OXALwOXmVlf0i0L\n99lHdlDM/Ba4K3BNg3w3UcyWfgb4prtWt8XATyj+OAa5B/CZBuUueS9J6ymWM9cOca+91PVViX/O\nlwE/CzzEzK6XdAzFknqYfrka+KKZ/coi133ZVwNXmdlRS9yzksfMfgA8TYUu+OvARyUdbGbbh6jv\nXtbqzGafMbPrgM8BfyXpAEkdST8j6VFLZPsI8FJJh0m6E/CHmTTvBc4A5szs3zPXoVhLTwM3AvMq\nROu9PwdLOlnSvcsZ0m1AD+hn7rORQi+4EZiQ9CrggCXqvxSXUszKupJOotAmajGzPnAW8OZSmOyq\nEIIXnq9PoRvkOBe4j6SnqxDZnwocDXx6iPp/EHiupGPKsv8P8DUz+9EQ9xpkyb5ahI0Us5NbJR0E\nvHo/yv80RRs9S9Jk+e8XJN23vP5Tqu37dWBrKVjPlv1xf0m/sFgBkp4p6c5lX95ans69b/vEWh1s\nFn5lWPj38Yb5nk3xMl1OoXV8FLjbEunfSTFAfYvim+pcij/23kCa91H8ivGPi92knLm8hGLwuoVi\nBnTOQJKjKMTnbRQaxtvM7AuZW32WYhbwfYolxC7SKXJTXkqhFSxM0xNblCX4A+DbwEUUIvmfA51y\n2foG4MvllP6hg5nM7GbgZIqZwM0Us7STzeymfa28mf0L8KcUs5DrKGZap+7rfTL3reurHH9NIbze\nBFzIcDO1wfJPpHiWn1AsKxfEaCi0qqPL9v2EmfUo2vQY4KqyDu+iWF4vxknAZZK2AX8DnLqwhNwf\nFn6FCZaB8lvu7WY2OH2fpRBPjy2np0Fwh2StzmxaoZyWPr6c8h9GMT32s6gXAhfFQBPc0RnpYCPp\npNLo6IrcLyiSpiV9uLz+NVUNl15Znv+epMeOsp77gSh+IryFYhn1XQp7huKi9COK5cjLxlG5IFhJ\njGwZVf6s933gVyh+xbgIeJqZXT6Q5kXAA8zsBZJOBX7NzJ4q6WgKge94ip8w/wW4T7n+DIJgFTLK\nmc3xFAZ0V5rZHuBDwCkuzSkUtihQiLG/XP7icgrwITPbbWZXUdiYHD/CugZBMGJGaWdzGNVfQq6h\nsGbMpjGzeUm3UdhCHEah2g/mzRpjqfAHOR2gMzH14JkDD8klC4JgGdi9bQvzu7YPYx+0+o36zOxM\n4EyA9ZsPt6NP/r0x1ygI1i6Xf/otQ+cd5TLqWqoWtXcntd7cm0aFa/uBFPYVTfIGQbCKGOVgcxFw\nlKR7SpqiMELyxk/nAM8pPz8F+LwVivU5wKnlr1X3pDBq+/oI6xoEwYgZ2TKq1GBeTGHV2gXOMrPL\nJL0WuNjMzqGwdnyfpCsoLE1PLfNeJukjFJa888DvDPNLVGcu/aWtu8eWPNZ8mkd9l6bX4Be8TnVZ\n2++my1ybcGkmq8e9yTRPf7J63Jvy90ir0nflWLd63XJfOb7ozCPLGbCr54/TTJ256nHXHXcy7e/7\nseM813J5fJ81wXyfTWTa3/3FDNNn6XF9OUkfNZkm5PrM9ZFvS2W8Agf/RmwotaZgpJqNmZ1LYcI/\neO5VA593Ab+xSN43UJi2B0GwBggL4iAIWiEGmyAIWiEGmyAIWmHV29ksRTcjEE9ur6qaE9urilh3\nV6qQaXdVxdR8/dYe1qmO4zbdTdL0Z6pK4fxsNY02pHnmuu6+7uvCC8YAvRlXrhcoc29BE4G4Rmzs\n7k7rIu8e4wXiTJ9N7DR3XC24syftD99HSbkZ+hPVxuxPpe3fm62mmVvn8jQQlednnajs+geg57Zg\n709W658V9R3KPLIX6DuujyZ2kVC5z34IxDGzCYKgFWKwCYKgFWKwCYKgFda0ZpMz9prYWRUWJrZW\ngx1oW7r7obZXz9ncXJImKXvCNe1MGgets362eqJfXbzbRPpdML9kODWwTI/2ptzxjFv/ZwwBE00g\nt/53Go259X9OM+jPLS0GdTKmm91d1ZPdHU5n27EnyaM9rnK9jM6mal06k9XG66/LNAzVxvQGefMZ\n/cWcQWff9cf8bNpQPXeuP+36bKJeg1Ivo5ntqZ7ruvdFGau9zmCe0GyCIFjpxGATBEErxGATBEEr\nrGnNJudUpjlnf7Gzut7Xth14+lu3VY5tT0azMXdfp9lo/To88prBVDVPZ13aPZ2e/35wTpaZr4/E\n8TJx5kzX/z5PTn/xjos4W5ZUn6m3D8k5UHrtreNsoXwfAmi3OzefEYNc+2vKNUymrt3JasP05qqJ\n8g6gro98s2X+CvuuT2zKaU6TOcMnlyej2Xj7L/Wqz5OzuRrss/3ZRDhmNkEQtEIMNkEQtEIMNkEQ\ntEIMNkEQtMKaFohzJAJevyq8Wc74q1cVF20uFSRry53LqNVOtJSryzC7zeULrx4mdlu5r5yOExtz\nemSiVbsbD2MAltU9bclj34dA0rY2n3Gw9fV1Tq5Zh9vk/fF1TbPUqqq5dqrpI3UaCMS5bfWSHf/c\nfXybLFa/IYiZTRAErTDu8Lu/L+lySd+S9K+Sjhi41pN0afnPb5QeBMEqY2TLqDL87lsZCL8r6ZzB\n8LsU8bGPM7Mdkl4IvAl4anltp5kdM6r6BUHQLqPUbPaG3wWQtBB+d+9gY2ZfGEh/IfDMEdanKNMb\nok24Dasm0g2TbLJq7KWcruON+qact503GMuUnWy45evahAaagdcV/CZY5dn6+/ZdGt8sw0hOmUc2\neaM4Z4zXyUzQvfFagzTJ8TD6RSPNydU/Z2/oDPK8gZ5vgyKTO5dzxOz5shv08zJJh6NcRuXC72ZD\n6JY8Dzhv4HhG0sWSLpT0pFFUMAiC9lgRv0ZJeiZwHPCogdNHmNm1ku4FfF7St83sh5m8e2N9T63f\n1Ep9gyDYd8YdfhdJjwH+GHiime3dXMbMri3/vxK4AHhQrhAzO9PMjjOz4yZm1i9f7YMgWFbGGn5X\n0oOAd1AMNDcMnN8kabr8vBl4BANaTxAEq49xh9/9C2AD8E+lgdV/mdkTgfsC75DUpxgQ3+h+xWpW\nh5wH66SPeuDE39l0q7Vk07rpqSRNYvXmvb4zO/XZbPVcf8btFDeZfhfUicY5o7IkLK5zWldG9UxE\ny0bRFZwXe074rAlMkXu+JHywi1RhvdyOeq7c+fRlSMryO/XNZPJ0vajsC8qUXROmuJPxju9UN5DE\nvAVlzqPbdX5up76O26kvibaQ2ylhoP77Y9837vC7j1kk31eAnx9l3YIgaJewIA6CoBVisAmCoBVW\nxE/fo6LfTVeYPbfeV7KDfmZHPafrZHd985pN1+2ANp02tY+I2Z9xeSYz+kVqc1ghZyCWrMPduj2n\no1jHG6LlynLrf+efmt0psca5NLeTn9euehktJSnH7ainXr1lmrk82YiY0y4C5qQ3tsvduHqYRqVM\ns3S94aLrIx+xIVd2rl+TqKW7ltZwwL1T+2HgFzObIAhaIQabIAhaIQabIAhaYU1rNrl1rV9zQ1U3\n8XY4ANrTwBEzKdvZ82SiW/anvBbhjqcz6/IhNJuu11K8I2ZGW2liUeE1Aa8H1K7/M+R0tr5rh3mc\nnU0mT2feNVR9l2HOnidn59RzfeZtgJrY2SS6SYO92LwtTt17kCsXcrZR7jhTl04DvasJMbMJgqAV\nYrAJgqAVYrAJgqAVYrAJgqAV1rRAnDOA86Kr3x0vZ0iXOiU2EMz87nKZunhx0ZediI/kQt5WD32o\n2rIy1SMvEua+cnzROaM+72DojzN1qROIc+3USwzn3HGmndR3D9Wgy3zbZsPidpfuo9wOeolAPOcq\nk7EE9KJ9EkK5wTQha4iZ9JGrW12fhVFfEAQrnRhsgiBohRhsgiBohTWt2WTX/06M8GvwnKNgsvZt\nsm6ti0JJRiPwMkOm/rXOdjlHTK8xeT1gyB2REuNAX5ecUVnt5lnpucTZsbP09SJRvSNpksX3WWYj\nr6SPhtk8yxv1Zd45b2yX9NGwu1jV9FmdIWCTdlyMmNkEQdAKMdgEQdAKMdgEQdAKI9VsJJ0E/A3F\nhufvMrM3uuunUWx6vhDi5Qwze1d57TnAn5TnX29m79nX8rNrbq95JDYz+7Ol876RjWo4SIOqJGvu\nzJp6WE1mX2m0nq9J06jPOkOIZstEbZ9lM1UPk43gc3uxraI+a8q4Y30DfNjMXuzyHgS8miJwnQGX\nlHlvGVV9gyAYLaNcRu2N9W1me4CFWN9NeCxwvpltKQeY84GTRlTPIAhaYCXE+n6ypG9J+qikhQia\njeOESzq9jAl+8fyu7ctR7yAIRsC4BeJPAUea2QMoZi/7rstE+N0gWBWMUiCujfVtZjcPHL4LeNNA\n3hNc3gv2tQJZkS0xtqtX4oYR64Yyfmri31mXJmeUldxjeRS/ocTSmizZtk58KkfTZ57lEk9rjUJz\nDpPJPVZunzVl3LG+7zZw+ETgu+XnzwInljG/NwEnlueCIFiljDvW90skPZHCgH4LcFqZd4uk11EM\nWACvNbMto6prEASjR7ZM07OVwPrNh9vRJ//e3uNGU+kGaVbVMqrRPVbulLzJ0rcJK3oZNdQ9Vkaf\nXf7pt7D9pquHat017YiZe0mHcWgb5sVNTAVz74p3iqu7SYbEcTS3/l+Glz3rYOgyNmnbRG9pkqcm\nzXINUEm5TdrS92HOqLLuHhmSwaVBOY2oy9ikz4Zk3L9GBUFwByEGmyAIWiEGmyAIWmFtazZDMLL1\nf4PbNnIJTfbKdpt/5UpaDm2xic7QJMswbVmn0WS+Mpflx4EGOlvSH010Hm8LlSunv3S/rsafdWJm\nEwRBK8RgEwRBK8RgEwRBK8RgEwRBK4RAPAxNDM+aMIzY6L4e/E59louuWLc74TCGgBmWLQJAHUMY\n9eWiNgy1U2KdTVym05KIBXVRKbJ5fBSQXNlLlzM0q8ARMwiCYC8x2ARB0Aox2ARB0AprW7NpsClR\nI+fBBlEPk0iJiVCSqV8dWS2lxkGvwfpf3mAsF6WywX1r6zYqGvSZj3SajS7q+zXpwwZlJzfNaGZ1\nzps1USiLNLbk9dx9mmlB7voIuzBmNkEQtEIMNkEQtEIMNkEQtMLa1mwaMIxmYxNpor7TBLxGkLPz\naKQR+Or5RH7N3WAt35n3eTK2IV4zyNzXR3b0Nj7Ltf5PHtkde30GwNyb3c+86Umamj6EBvpdA9uc\nRDfJaCu+bdWr3tj3IUBnzudJOyDt+/q6LJe9zkhnNpJOkvQ9SVdIekXm+lskXVr++76kWweu9Qau\nnePzBkGwuhhr+F0z+72B9L8LPGjgFjvN7JhR1S8IgnZZSeF3nwZ8cIT1CYJgjKyE8LtIOgK4J/D5\ngdMzZVjdCyU9abFCIvxuEKwOVopAfCrwUTMblKuOMLNrJd0L+Lykb5vZD31GMzsTOBOKUC51BdUJ\nwl4kBOhPVRP1J9M0vSlXTiOB0u2+1sB4MK2cy5ITcueduOiExM5cWlBnj0+TKXuuxtAsJ1YPIzbW\nOFX6toa0P/pTaRp/rjdVrVyun2v7rFO/U2JiVJlrp7ml+6y7p77Pcmks6ddq5XLCc65+wzDKmU1t\n+N0BTsUtoczs2vL/KylC7z4ozRYEwWphrOF3AST9HLAJ+OrAuU2SpsvPm4FHAJf7vEEQrB7GHX4X\nikHoQ1YNzXlf4B2S+hQD4hsHf8UKgmD1MVLNxszOBc51517ljl+TyfcV4OdHUim//ncaTX8yXef2\nppc+BujNmkvj1v/T6VreJt25CW99l5aTGFi59T+9jCOgW7t3d1Vv3NmdFjOx0+XJpPG7cnX31G/K\nlWhMSbSIXDkuTbeBhub7bF1amd6M67NZV7np1MKtM10VMCa61TTKaTaOfq/jjjPv3O7qiyl33NuZ\n5um6c7Yjk6bGEzmJsEozJ9wmhLtCEAStEINNEAStEINNEAStEINNEAStsFKM+lqjzmM4Z/zlxcb5\nDamINr++qqLZuqqQOLkutYqbnqmem5qoWlRNduuVuTknNu6eS9XS3buq5+Z3Vrtd21NLRpuo3tc6\nOXf46mHiXd5kp7gGJH3m3lpvwAepYD+/Pu2z/gbX3huqFm8b1qWq+IZpl2aqmmbKu2tn2NWrPsCO\nufQBtu6qvnTbd1SP57an/Wxd1zBZo9ClvcctY9Q3VOSQDDGzCYKgFRoNNpIeKukiSdsk7Sm3f7h9\n1JULgmDt0HRmcwaFV/YPgFng+RTbRwRBEDSisWZjZldI6pbOkv8g6T+AV46uastAbq3pnficXNFs\n/Z8RHg6s6i8bNu6qHB+8fkeSZfPstuotJqt5Zrs578cqO3vVtfttczNJmpt2bqgc37x9XeV4+0Sa\nZ17V+6qffi95B89+4jzYJLxFmsSTOF4mhphpHt9nXp8BmDmwqrfc+YBqfxy64bYkz11nqhP6zZPV\nPOu8N2SGHU4YvGluQ5Lm+l0HVI5/su3AyvGNE2meXe6410v/vDvOwNPrX50RCitNB5sdpX/TpZLe\nBFxH6D1BEOwDTQeMZ1H4N70Y2E7hzf3kUVUqCIK1R6OZjZn9uPy4E/iz0VUnCIK1SqPBRtLJwOuA\nI8o8AszMDlgy4wok2anf76if24jJOVF6GxqA9Ruq6/9DD6iu7e+x4ZYkzxGzN1eOD5ms5lmf85B0\nbO9X7S9umEu75MdTB1eOp7ubKsfXZgwpts9XJ7293Rmnvl1LR53MRpSoo5EjZvW4P5WxoXFOlhPr\nUs3moI3VnR2PPKDaH/dZf0OS557TN1aO7zpxa+V4Y8crJylb+1WN7Pr5OyVprpq8c+V4yhnE9DJ9\ndsOcc9bcldpP9Xe5d9l7Zta1/37Y3DTVbP4a+HXg224riCAIgkY0/e65GvhODDRBEAxL05nNy4Fz\nJX0R2Du3N7M3j6RWQRCsOZoONm8AtgEzQMYSJQiCYGmaDjaHmtn9R1qTtqjZqT8XBcELkH63NoCN\ns1Vh8M7OYM+LwQD3nv5p5fjQyaqIfIDqBeLbrSoQH9DZWZvHGwLevifdenDnjup3Sn8qIza6thoq\nOkQDTDVCdC5Mrtv1cGo6NZDcNFNtq7s5gz0vBgMcNXV95fhwt4Xhxk79n9TWftVYsJGo3KuKylt2\nr0/S3Do9WzneMZlz1nQOtr4tleu05VFPmmo250o6cV9vLuksSTdI+s4i1yXpb8vwvN+SdOzAtedI\n+kH57zn7WnYQBCuLpoPNC4HPSNop6XZJWxs6Yp4NnLTE9ccBR5X/Tgf+HkDSQcCrgYdQRNZ8taRN\ni90kCIKVT6PBxsw2mlnHzGbN7IDyuNbGxsy+BGxZIskpwHut4ELgTpLuBjwWON/MtpjZLcD5LD1o\nBUGwwmnsiCnpAcCRg3nM7J/3s/zFQvTuS+je0ylmRUytH2Lyk0TIzKxP3ZDcyeygP+l2ivJOlDkH\nvY3dqmbgNZqN2TCUDucTurWbaja+bF83X3dIn7GX+1rybZVd7y8DNbfNbu7k6tvNbETmDeV8O+WM\nKjeq2nZeo9nQSZ1aU6oazcZ+2s++bF83X3fIPGOmz2o3whpRF0JzC+KzgAcAl/Hfr7cB+zvY7Df7\nGn43CILx0HRm81AzO3oE5S8Wovda4AR3/oIRlB8EQUs0FYi/KmkUg805wLPLX6UeCtxmZtdRRNE8\nsQzDuwk4sTwXBMEqpenM5r0UA871FBbEC46YD1gqk6QPUsxQNku6huIXpkmKzG+niJb5eOAKYAfw\n3PLaFkmvo4gXDvBaM1tKaB4ev/DKLWrdUrjvo1ACc25zKW/L4jdMAtjaq9pF3O7tLRpsDu7tbPw9\nc2X7uvm6Q+YZc3VJvFoXr+d+UXPfnMzm69vLiE57nKGQbyfv5Aqw1aptt7XvdZ0mjphVvWWrpeX4\nsn3dfN0h84y5Defr+miEQkTTwebdFHvafJtGfwIFZva0musG/M4i184CzmpaVhAEK5umg82NZnbO\nSGsSBMGapulg8x+SPgB8iqoj5th/jQqCYHXQdLCZpRhkBl0WVsRP30EQrA6abgv63FFXpDXqojhm\nIgL6Hel7u1PPv607q8ZcN05Vd7+f7VZ3y8txe78q7i7bTn07q2Xf6KIt+LoD9N0zdvekorhvqyTa\n5TKJjXLbKKknd5zJ43Ya3LM7dUq8ZVe1va+brrbdxm51t7wcW4faqa/qRJndqW93tezrXLQFX3dI\nn1FzqSju2yo5HuGWVU2N+maA5wH3o9hmAgAz+80R1SsIgjVGUzub9wF3pfBZ+iKFkd3WUVUqCIK1\nR9PB5t5m9qfAdjN7D/CrFB7ZQRAEjWgqEC94it0q6f7A9cAho6nSaEl8B/2aNaPZdF1kgf6OVLPZ\nMVnVTn6i6hp7dyY64S17quvucUXE3LE1NSqTe0bfBgDeTzRpy8YWWQPkJIOaPvOaGkDHRX6Y35G2\n/5bJqnbSdS9HznDOR68cV0TMLVvTzbP8M/roF5C2VaJ31bX/fkg6TQebM0u3gT+hcDHYAPzp8MUG\nQXBHo+lgcyClKwHw1vL/eUnHmNmly1+tIAjWGk01mwcDL6DYU+ZQiv1jTgLeKenlI6pbEARriKYz\nm7sDx5rZNgBJrwb+H/BLwCXAm0ZTvf2kgYOeX7N2M0tu2+k33E7H6Hmq2sk2Z+Owe1dq53HzTFU7\nmZqoCkaTmQ2fPHPO+W73XFqOL3t+Z7XbtT3VoCa2V+/b3Zlb/7vjedfgDfSXJnjtJ9FsMtKWr691\n01fdW8Rc7+x3tu1OnWd/Ml3VTjZMVW2hpjoZox9frtPvdsxlnHR3VXW07Tuqx3Pb037ubG2g2bi2\nqrWVWkaazmwOYcBNgUIwvouZ7XTngyAIsjSd2bwf+JqkT5bHTwA+IGk9cPlIahYEwZqiqbvC6ySd\nBzyiPPUCM7u4/PyMkdQsCII1ReMNz8vB5eLahEEQBBkaDzZrBW/U1+lVT1jGQCzVg9M0ckJtb7c3\nKktF2DlnCMhE/e74icDqd9TrZermnmliV/XGOX/PCSewdjNp/LlRiY1Jn7lysqJ+x7VDZgPGnhNq\n5/ZU2+WW7alwe5uLhuojGijXe2AqAAARcklEQVQTecPTd+9KP9Nn3hFW3jE2I9j7cxM7Mmlcn3kj\n1lyf1e7u15CmAnEQBMF+MbLBpkHo3WeUIXe/Lekrkh44cO1H5flLJcXSLQjWAKOc2ZzN0lEsrwIe\nZWY/D7yOMvbTAI82s2PM7LgR1S8IghYZmWZjZl+SdOQS178ycHghheHg6Klz6pvLLFBdpMfsutat\nu71BVcanD5twelHH6TpNohPWGLwBdOarN0oMuzJGcYnBXs5wbs5vauXqMsT6P3fd7+ckp7N15jLa\nhJdsMlEzfL6Od7idTHU2m6h25Lz/us5pNsmGbfWbf034urn2z25mtsenSe9bZ4iZ1dnWmGbzPOC8\ngWMDPifpkjK87qJIOl3SxZIunt+1faSVDIJgeMb+a5SkR1MMNo8cOP1IM7tW0iHA+ZL+08y+lMsf\n4XeDYHUw1pmNpAcA7wJOMbObF86b2bXl/zcAHweOH08NgyBYLsY22Ei6B0V0hmeZ2fcHzq+XtHHh\nM0VEh+wvWkEQrB5GtoxqEHr3VcDBwNtUCLDz5S9PdwE+Xp6bAD5gZp9Ztnp5x+TEyC+TaY8X0TKG\nc844ypy2aLlh3Rme5SL/JuUk4YLd9VykgZoIEl5wzd0nKzzX7PK2XMZgdRExsrsruuOc8OnbwQvE\nvg+hgbFgrg9r6p+LMevbNvkhIxcFJNk5Me2AtO9r6raMjPLXqLrQu88Hnp85fyXwwDRHEASrmZXy\na1QQBGucGGyCIGiFsf/0PW7qNBzIaAQZwz+/9jXVOwLWkquLr2CdHpA5V3fc9L61dcvddoh2qdOp\nvDMtQN/dOGtr53UQ99Wb1dDq6pvts6XT5I1EfZqlDShz92nyLoxMZ8sQM5sgCFohBpsgCFohBpsg\nCFphbWs2mfV1rS1LA82mkS1FA++1urV8k/V/3Ro8l6ZO98nWrQFJ2w6jUzUqqOaYzKZoOc3G60ee\nBv2cZMkUVKuTDKOtNHFyXS79ZZn6MWY2QRC0Qgw2QRC0Qgw2QRC0Qgw2QRC0wtoWiEfFMILqMMZe\nTcTqJE8DgbJB3YbB64jZ2y6H2NignfyprBFcXacN0WfjFPWXhdyPKsukEMfMJgiCVojBJgiCVojB\nJgiCVgjNxtFk/b9c963VHhoUPM61vGdkRn2JwaS7nDNwG4E2BMPpbMP0c63D7ZD93GSDtlERM5sg\nCFohBpsgCFphnOF3T5B0Wxli91JJrxq4dpKk70m6QtIrRlXHIAjaY5SazdnAGcB7l0jzb2Z28uAJ\nSV3grcCvANcAF0k6x8wu3+caNFlzN1nDDrHOHWpNPYzzo9+Au4nNSYNNrhqVPYwwMowjbE2arL/k\nMmgTjfqwka6272X7tvWOvcM+30j6rCEjm9mUQeW2DJH1eOAKM7vSzPYAHwJOWdbKBUHQOuPWbB4m\n6ZuSzpN0v/LcYcDVA2muKc9lifC7QbA6GOdg8w3gCDN7IPB3wCeGuYmZnWlmx5nZcRMz65e1gkEQ\nLB9jG2zM7HYz21Z+PheYlLQZuBY4fCDp3ctzQRCsYsZm1CfprsBPzcwkHU8x8N0M3AocJemeFIPM\nqcDThyqjkYHVCLeTr2GYneJqhcHM10eaZzSWXcMIqomTYu4edU6JLbIsfTaMgeQK6rNhGWf43acA\nL5Q0D+wETjUzA+YlvRj4LEUU1bPM7LJR1TMIgnYYZ/jdMyh+Gs9dOxc4dxT1CoJgPIz716ggCO4g\nrGlHTB9FEHKRBhvkGWbXem9rl1vL+wiM/ribydOt3qhRHv+V0qBuTajdFKpJdM6ayI+N8mS1OatP\n47P4dvEGk4C50Jpp+2ca053q+z7K6Wx10TmHlXBq+qxJdM5hiZlNEAStEINNEAStEINNEAStsLY1\nm8xasztXXbR2/PF8vc6TDa+YZHLaSkZL6U9U0/Qn3XFGTPFr+X4DzcaX49Mkmg402r28br2vTFt2\n5tyJGg2tyOP7yF3P9ZnXfppoNk6jsYk0U7+7dJ/1mugv7q/O909xzuVp0meOJhvmq0FbDha1Pw6u\nMbMJgqAVYrAJgqAVYrAJgqAVYrAJgqAV1rRA3MkJxLvNHffdcZpJe6pp1KsLMQnWrY7jNpGO6/2p\n6rneTPV4fjaTx/eYE6KzYuOUO550xxlRuQlebEyE24zTYiLcOsG4kajv+3BP2h+deXeuvsuwRLBP\n27/n+iwx2Mv8RXnhORGVp9M8PddniUDcoM+aGOj5PrPdNX0WAnEQBCudGGyCIGiFGGyCIGiFNa3Z\nqJcaKHmNZmJbVTTo7vBWZ6Dd7tx8Rljwhn7d6qLaptOm7s9UxRM5MSXnCNibXHrRnFvL+/V/ouFk\ndYZ6R0b1nBHcHn+TNI/XCJLrmT7zGs3Ezmr7d3elN9VcNU3uXfDYZLXx+lM5C0nXrxNeTGlgiOna\n2/cPQG/G5ZlyDqA5h1tXdNZA0jVVd1d9dI7+3PJs3BUzmyAIWiEGmyAIWiEGmyAIWmGUexCfBZwM\n3GBm989c/9/AMwbqcV/gzma2RdKPgK0UbnrzZnbcMHXIrf+9HY3XaDq370hvtHNX5dDmUl0n0Wwm\nqk3bmUmNKTQ/W72Fs83pTKcL85ztUPUe6blUI6jW1dvdFAX5G6dJUh2kurbP1TXVoZw2lHXEdLZQ\nTqPp7PBiEWhXtY+U0dl8XTTpGmpdKqZ0vc3MdLWhspt0eVscb+eUs7OZsSWPyTiJmivca2oAtmdp\n/SWbZ/CdWqF2NmcDJy120cz+wsyOMbNjgFcCXzSzwQiajy6vDzXQBEGwslgp4XefBnxwVHUJgmD8\njF2zkbSOYgb0sYHTBnxO0iWSTh9PzYIgWE5Wgp3NE4AvuyXUI83sWkmHAOdL+s9yppRQDkanA0yt\n3zT62gZBMBQrYbA5FbeEMrNry/9vkPRx4HggO9iY2ZnAmQDrNx9eUcj8LmQAcmJjYrDnxGCA/rbt\n1TJ37c5UxN13yomLvYyDZ8cLwlXlsLMu7R75rfmcYpeN4uDtzib9cUZs9EJzTiCedwKrE4z9rnZQ\nv8NcLrqC3z1OTuT3YjCAdro+ms8Y/rn2xznYdvx1oO+NM72gnWmnuggY/Uz796edEd+0K2gi52Xp\nysmIvX33TF4Qzht4DnxOLzdmrMsoSQcCjwI+OXBuvaSNC5+BE4HvjKeGQRAsF+MMvwvwa8DnzGxw\n6nAX4OMqtieYAD5gZp8ZVT2DIGiHsYXfLdOcTfET+eC5K4EHjqZWQRCMi5Wg2bRKogk4Yy/LOVk6\nIz6bS43IasudyDS1K0t9p/tk9Iv6gurPJZEfs5E366MTJAZ43mBvGAOwXDlJdEuXqJ/RL9w5y6RJ\n7+PFrZwAk56q3jSTxUfaGKb93bEyRn3eotCUqiTm3qmknMyGZ/tjyDfI2H/6DoLgjkEMNkEQtEIM\nNkEQtEIMNkEQtMIdTiBOvI69kVM3I6q5Xfc0mdlazRv1eUHYexQDuF3ezNUlt1PfUHgN0GuLucgD\nTijMejP3a3Z5G8YCrJHA6uqWMb7zbZttSS+G+t0VMxExEhHcG+w1EOgTcu1U00fm2z5XTq5f/bm6\nPlzs3BDEzCYIglaIwSYIglaIwSYIglZY05qN5ZzKXJRDm/UOk+uSPD6yY3anPkei2WR26rN11S30\n+zPVPLmIjDnnxkq5DaIgylW/mxEVEofJBpESOm4XuNxOfbn6VcvNOA/6SJWuneTDRQB47S0XxdT3\nq9PVfPQLgJ7bPdFHt2yis3n9K+cw7Nuy7wz0vHFetpyMI6Z8H7l3IRf9YrDP9kdFjJlNEAStEINN\nEAStEINNEAStsKY1G7/WB5ifXXrzo+5k6hWndVW9RfM1wgPpJkWWiZTgNYH52Wqa3kzNTvcZcuv/\nrvMblTMG6eckqKWDIBRJnCaTajhpHr8Rlsf7QgL0ZpzNjPcRzdjD+D5KnC5zZbv75CJi9maraZII\npTnzF7epmNdjJjI6j5z9S3+Pc6BsME3I2UYlGs3uajndzL5wlT7bD5ubmNkEQdAKMdgEQdAKMdgE\nQdAKMdgEQdAKa1ogTsQ7gPVOBJysirSaz0U0cAJxJqxvgg/Tmos04I3VXH1z9ffCoDeS6+5J6+bF\nxiTaQu4rp4lA7MpOjAcz7ZQTsAfxbQDgs/Sd8V1OdB5ml0NvkJf7gcFHH0iM+jKvnBfOJ3Y6wXg+\nzdR1QT5SB+K0nIQhRP2sgeGgqLwSBWJJh0v6gqTLJV0m6aWZNJL0t5KukPQtSccOXHuOpB+U/54z\nqnoGQdAOo5zZzAMvM7NvlKFZLpF0vpldPpDmccBR5b+HAH8PPETSQRTRGI6jGEsvkXSOmd0ywvoG\nQTBCRhnr+zoz+0b5eSvwXeAwl+wU4L1WcCFwJ0l3Ax4LnG9mW8oB5nyKEL1BEKxSWtFsJB0JPAj4\nmrt0GHD1wPE15bnFzufuvTf8LrD74ve8bC0GtNsM3DTuSoyAtfpcsHaf7WeHzTjywUbSBuBjwP8y\ns9uX+/6D4XclXWxmxy13GeMmnmv1sVafTdLFw+Yd6U/fkiYpBpr3m9k/Z5JcCxw+cHz38txi54Mg\nWKWM8tcoAe8Gvmtmb14k2TnAs8tfpR4K3GZm1wGfBU6UtEnSJop4358dVV2DIBg9o1xGPQJ4FvBt\nSZeW5/4IuAfsjfd9LvB44ApgB/Dc8toWSa8DLirzvdbMtjQo88zlq/6KIp5r9bFWn23o55I18IYN\ngiDYX8JdIQiCVojBJgiCVlh1g42kkyR9r3RxeEXm+rSkD5fXv1ba+KwKGjzbaZJulHRp+e/546jn\nviDpLEk3SMraPy3lsrLSafBsJ0i6baC/XtV2HYdhf12NFsXMVs0/oAv8ELgXMAV8EzjapXkR8Pby\n86nAh8dd72V8ttOAM8Zd1318rl8CjgW+s8j1xwPnUbh+PhT42rjrvIzPdgLw6XHXc4jnuhtwbPl5\nI/D9zLu4z/222mY2xwNXmNmVZrYH+BCFy8MgpwDvKT9/FPjl8mf4lU6TZ1t1mNmXgKV+SVzMZWXF\n0+DZViW2f65Gi7LaBpsmbgx705jZPHAbcHArtds/mrpoPLmctn5U0uGZ66uNxq4pq5SHSfqmpPMk\n3W/cldlXhnA1WpTVNtjc0fkUcKSZPYDCOfU9NemD8fIN4AgzeyDwd8AnxlyffWK5XY1W22DTxI1h\nbxpJE8CBwM2t1G7/qH02M7vZzBb2v38X8OCW6jZK1qxripndbmbbys/nApOSNo+5Wo3YD1ejRVlt\ng81FwFGS7ilpikIAPselOQdY2GzrKcDnrVS0Vji1z+bWxE+kWEuvdhZzWVn1SLrrgl4o6XiKv7cV\n/8W3n65Gi7KqtgU1s3lJL6bwk+oCZ5nZZZJeC1xsZudQNNL7JF1BId6dOr4aN6fhs71E0hMpNibb\nQvHr1IpG0gcpfpXZLOkaik3RJmFpl5XVQINnewrwQknzwE7g1FXyxTe0q9FShLtCEAStsNqWUUEQ\nrFJisAmCoBVisAmCoBVisAmCoBVisAmCoBVisAmWHUlHLuYJvUj60yQdOso6BeMnBptgJXAaEIPN\nGicGm2BUTEh6v6Tvlk6j6yQ9WNIXJV0i6bOS7ibpKRSRT99f7vkyK+lVki6S9B1JZ64Sr/2ghjDq\nC5ad0lP4KuCRZvZlSWdRuFb8GnCKmd0o6anAY83sNyVdAPyBmV1c5j/Iyg3uJb0P+IiZfWoMjxIs\nI6vKXSFYVVxtZl8uP/8jhbn7/YHzy4lKF1jMl+bRkl4OrAMOAi6j8HgPVjEx2ASjwk+ZtwKXmdnD\nlsokaQZ4G3CcmV0t6TXAzGiqGLRJaDbBqLiHpIWB5enAhcCdF85JmhzYTGorxfaT8N8Dy03lfipP\naavCwWiJwSYYFd8DfkfSd4FNFJtHPQX4c0nfBC4FHl6mPRt4e+lhvBt4J/AdCg/4iwjWBCEQB0HQ\nCjGzCYKgFWKwCYKgFWKwCYKgFWKwCYKgFWKwCYKgFWKwCYKgFWKwCYKgFf4/lIfuHHCxXooAAAAA\nSUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1pWhCblZASbb",
        "colab_type": "text"
      },
      "source": [
        "By inspection we can see that the energy function has a number of interesting properties. First, note that the function is periodic in $\\beta$ and $\\gamma$ with shorter periods than one might naively expect given the definition of the gates. The details of why that's true will take us away from the main content of this Colab, but it's a good thing to understand so that the parameter space can be efficiently truncated.\n",
        "\n",
        "The other main thing to notice is that there are many local minima and maxima. This makes it challenging to use gradient-based methods for optimization. We'll see that explicitly next. Part of the challenge for algorithms of this type is finding efficient ways to optimize the paramters."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9MwLJ0LCC0un",
        "colab_type": "text"
      },
      "source": [
        "#### Gradient Descent\n",
        "\n",
        "For practice let's try to minimize the expectation value of the energy using gradient descent. We know that there are local minima that we might get stuck in, depending on initialization, but it's still a worthwhile exercise.\n",
        "\n",
        "The first step is to define a function which approximates the gradient of the energy. We'll do this by symmetric difference, i.e., $f'(x) \\approx (f(x+\\epsilon)-f(x-\\epsilon))/(2\\epsilon)$. You should experiment with different values of $\\epsilon$ as well as different formulas for the gradient. "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "_Ge8sxrab6hV",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def gradient_energy(gamma, beta, qaoa, h):\n",
        "  \"\"\"Uses a symmetric difference to calulate the gradient.\"\"\"\n",
        "  eps = 10**-3 # Try different values of the discretization parameter\n",
        "  \n",
        "  # Gamma-component of the gradient\n",
        "  grad_g = energy_from_params(gamma + eps, beta, qaoa, h)\n",
        "  grad_g -= energy_from_params(gamma - eps, beta, qaoa, h)\n",
        "  grad_g /= 2*eps\n",
        "  \n",
        "  # Beta-compoonent of the gradient\n",
        "  grad_b = energy_from_params(gamma, beta + eps, qaoa, h)\n",
        "  grad_b -= energy_from_params(gamma, beta - eps, qaoa, h)\n",
        "  grad_b /= 2*eps  \n",
        "  \n",
        "  return grad_g, grad_b"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xaRHKo2LDwdm",
        "colab_type": "text"
      },
      "source": [
        "Now we'll implement a gradient descent algorithm that minimizes the energy. Note that it will get stuck in local minima depending on the initialization."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "BYcCDiYDB8ef",
        "colab_type": "code",
        "outputId": "883215a4-8663-4933-9322-5739288e9459",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 161
        }
      },
      "source": [
        "gamma, beta = 0.2, 0.7 # Try different initializations\n",
        "eta = 10**-2 # Try adjusting the learning rate.\n",
        "for i in range(151):\n",
        "  grad_g, grad_b = gradient_energy(gamma, beta, qaoa, h)\n",
        "  gamma -= eta*grad_g\n",
        "  beta -= eta*grad_b\n",
        "  if not i%25:\n",
        "    print('Step: {} Energy: {}'.format(i, energy_from_params(gamma, beta, qaoa, h)))\n",
        "print('Learned gamma: {} Learned beta: {}'.format(gamma, beta, qaoa, h))"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Step: 0 Energy: 0.3555179724197741\n",
            "Step: 25 Energy: -0.60656794301801\n",
            "Step: 50 Energy: -0.6068781974520839\n",
            "Step: 75 Energy: -0.6068778212318446\n",
            "Step: 100 Energy: -0.6068778804908308\n",
            "Step: 125 Energy: -0.6068782185651193\n",
            "Step: 150 Energy: -0.6068782066830509\n",
            "Learned gamma: 0.19753242844183583 Learned beta: 0.26843952111680774\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6UpLOKeso7uo",
        "colab_type": "text"
      },
      "source": [
        "### Results\n",
        "\n",
        "We've optimized our parameters. How well did we do?"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_Pv7GZX0EJdr",
        "colab_type": "text"
      },
      "source": [
        "For a $3\\times 3$ grid we have $9$ qubits and $12$ interacting nearest-neighbor pairs. If all of the qubits are in the $|0\\rangle$ state or all are in the $|1\\rangle$ state, then the energy-per-qubit is $-12/9 = -1.33$ at zero external magnetic field $h$, and will be close to that if the magnetic field is small. Notice that the QAOA algorithm we analyzed above is __not__ getting close to that ground state. Is this a problem?\n",
        "\n",
        "Well, not really. The QAOA algorithm still succeeds if we can  find the ground state after a small numbe of measurements. The QAOA prepares a certain state which is a linear combination of the ground state and many other states. When we measure the qubits, we find the ground-state configuration with some probability. If that probability is relatively large, then after a reasonably small number of measurements we'll locate the ground state.\n",
        "\n",
        "Practically speaking, this means we should measure the state prepared by the QAOA several times and record the lowest-energy state we find. The QAOA can be successful by biasing these measurements toward the ground state, even if they do not produce the ground state with $100\\%$ probability"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9epP9VB8GGgD",
        "colab_type": "text"
      },
      "source": [
        "Let's make a copy of our qaoa circuit for measurement purposes and attach a measurement gate to each qubit:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab_type": "code",
        "outputId": "41819d0e-0c9a-40e2-9b45-f4da116447e3",
        "id": "bIoIBl-aSjKE",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 379
        }
      },
      "source": [
        "measurement_circuit = qaoa.copy()\n",
        "measurement_circuit.append(cirq.measure(*[qubit for row in qubits for qubit in row],key='m'))\n",
        "measurement_circuit"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<pre style=\"overflow: auto; white-space: pre;\">                                                ┌──────────────────┐               ┌──────────────────┐\n",
              "(0, 0): ───H───ZZ─────ZZ─────Z^(0.5*g)───X^b────────────────────────────────────────────────────────────────────────────────────────────M('m')───\n",
              "               │      │                                                                                                                 │\n",
              "(0, 1): ───H───┼──────ZZ^g───ZZ──────────ZZ──────Z^(0.5*g)─────────────X^b──────────────────────────────────────────────────────────────M────────\n",
              "               │             │           │                                                                                              │\n",
              "(0, 2): ───H───┼─────────────┼───────────ZZ^g────ZZ────────────────────Z^(0.5*g)────X^b─────────────────────────────────────────────────M────────\n",
              "               │             │                   │                                                                                      │\n",
              "(1, 0): ───H───ZZ^g───ZZ─────┼───────────ZZ──────┼────────Z^(0.5*g)────X^b──────────────────────────────────────────────────────────────M────────\n",
              "                      │      │           │       │                                                                                      │\n",
              "(1, 1): ───H──────────┼──────ZZ^g────────ZZ^g────┼────────ZZ───────────ZZ───────────Z^(0.5*g)─────────────X^b───────────────────────────M────────\n",
              "                      │                          │        │            │                                                                │\n",
              "(1, 2): ───H──────────┼──────────────────────────ZZ^g─────┼────────────ZZ^g─────────ZZ────────────────────Z^(0.5*g)───X^b───────────────M────────\n",
              "                      │                                   │                         │                                                   │\n",
              "(2, 0): ───H──────────ZZ^g────────────────────────────────┼────────────ZZ───────────┼────────Z^(0.5*g)────X^b───────────────────────────M────────\n",
              "                                                          │            │            │                                                   │\n",
              "(2, 1): ───H──────────────────────────────────────────────ZZ^g─────────ZZ^g─────────┼─────────────────────ZZ──────────Z^(0.5*g)───X^b───M────────\n",
              "                                                                                    │                     │                             │\n",
              "(2, 2): ───H────────────────────────────────────────────────────────────────────────ZZ^g──────────────────ZZ^g────────Z^(0.5*g)───X^b───M────────\n",
              "                                                └──────────────────┘               └──────────────────┘</pre>"
            ],
            "text/plain": [
              "                                                ┌──────────────────┐               ┌──────────────────┐\n",
              "(0, 0): ───H───ZZ─────ZZ─────Z^(0.5*g)───X^b────────────────────────────────────────────────────────────────────────────────────────────M('m')───\n",
              "               │      │                                                                                                                 │\n",
              "(0, 1): ───H───┼──────ZZ^g───ZZ──────────ZZ──────Z^(0.5*g)─────────────X^b──────────────────────────────────────────────────────────────M────────\n",
              "               │             │           │                                                                                              │\n",
              "(0, 2): ───H───┼─────────────┼───────────ZZ^g────ZZ────────────────────Z^(0.5*g)────X^b─────────────────────────────────────────────────M────────\n",
              "               │             │                   │                                                                                      │\n",
              "(1, 0): ───H───ZZ^g───ZZ─────┼───────────ZZ──────┼────────Z^(0.5*g)────X^b──────────────────────────────────────────────────────────────M────────\n",
              "                      │      │           │       │                                                                                      │\n",
              "(1, 1): ───H──────────┼──────ZZ^g────────ZZ^g────┼────────ZZ───────────ZZ───────────Z^(0.5*g)─────────────X^b───────────────────────────M────────\n",
              "                      │                          │        │            │                                                                │\n",
              "(1, 2): ───H──────────┼──────────────────────────ZZ^g─────┼────────────ZZ^g─────────ZZ────────────────────Z^(0.5*g)───X^b───────────────M────────\n",
              "                      │                                   │                         │                                                   │\n",
              "(2, 0): ───H──────────ZZ^g────────────────────────────────┼────────────ZZ───────────┼────────Z^(0.5*g)────X^b───────────────────────────M────────\n",
              "                                                          │            │            │                                                   │\n",
              "(2, 1): ───H──────────────────────────────────────────────ZZ^g─────────ZZ^g─────────┼─────────────────────ZZ──────────Z^(0.5*g)───X^b───M────────\n",
              "                                                                                    │                     │                             │\n",
              "(2, 2): ───H────────────────────────────────────────────────────────────────────────ZZ^g──────────────────ZZ^g────────Z^(0.5*g)───X^b───M────────\n",
              "                                                └──────────────────┘               └──────────────────┘"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 14
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_OOqzrQwGTJZ",
        "colab_type": "text"
      },
      "source": [
        "Now we'll instantiate a simulator and measure the output of the circuit repeatedly:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "KbIu8eyNSK_t",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "num_reps = 10**3 # Try different numbers of repetitions\n",
        "gamma, beta = 0.2,0.25 # Try different values of the parameters\n",
        "simulator = cirq.Simulator()\n",
        "params = cirq.ParamResolver({'g':gamma, 'b':beta})\n",
        "result = simulator.run(measurement_circuit, param_resolver = params, repetitions=num_reps)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "EudMLjzNGadh",
        "colab_type": "text"
      },
      "source": [
        "Finally, we'll compute the energy for each of our measurement outcoems and look at the statistics. We start with a helper function which calculates the energy given a set of measurement outcomes:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Oa6kAObJTZRi",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def compute_energy(meas, h):\n",
        "  Z_vals = 1-2*meas.reshape(n_rows,n_cols)\n",
        "  energy = 0\n",
        "  for i in range(n_rows):\n",
        "    for j in range(n_cols):\n",
        "      if i < n_rows-1:\n",
        "        energy -= Z_vals[i, j]*Z_vals[i+1, j]\n",
        "      if j < n_cols-1:\n",
        "        energy -= Z_vals[i, j]*Z_vals[i, j+1]\n",
        "      energy -= h[i,j]*Z_vals[i,j]\n",
        "  return energy/(n_rows*n_cols)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kkUl5LYnG7E7",
        "colab_type": "text"
      },
      "source": [
        "Now consider the 10 most common outputs of our measurements, and compute the energies of those:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "t2SHZj_-TTFS",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "hist = result.histogram(key='m')\n",
        "num = 10\n",
        "probs = [v/result.repetitions for _,v in hist.most_common(num)]\n",
        "configs = [c for c,_ in hist.most_common(num)]"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "6-jbvrc_WOgP",
        "colab_type": "code",
        "outputId": "9bb8140c-d84c-43d2-efc0-ca86e76dc4bc",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 563
        }
      },
      "source": [
        "plt.title('Probability of {} Most Common Outputs'.format(num))\n",
        "plt.bar([x for x in range(len(probs))],probs)\n",
        "plt.show()\n",
        "meas = [[int(s) for s in ''.join([str(b) for b in bin(k)[2:]]).zfill(n_rows*n_cols)] for k in configs]\n",
        "costs = [compute_energy(np.array(m), h) for m in meas]\n",
        "plt.title('Energy of {} Most Common Outputs'.format(num))\n",
        "plt.bar([x for x in range(len(costs))],costs)\n",
        "plt.show()\n",
        "print('Fraction of outputs displayed: {}'.format(np.sum(probs).round(2)))"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAEICAYAAACzliQjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAFohJREFUeJzt3X+0XWV95/H3pwmCggWFOEoSCTZo\nJ/ijagSsrXWZtgNjNbWFTtAqOlScaVE6Oto4a6pI/dl2RDuia1FRkaiAUZxUorgc2unSodGAVo2I\nRuRHwB/hhygIhuh3/tj7DsfjvbnnJvfeE/K8X2tlcfazn3P2d+9z7uc859n7HFJVSJLa8EvjLkCS\nNH8MfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6e7EkZyZZt5v3fVGSz+5i/SeTnDJZ3yR3JnnU7mx3\nhjU+MMk/JLkjyUfmenuSDP1Zl+S6JHf3wfm9JO9PctC46xpWVSdU1flTrDuoqq4F6Ot/wxyVcSLw\nb4BDq+qk4ZVJHpvksiS3JPmFL5QkeWiSS5LcleT6JM+bakP9G2glOWOo/Yy+/cw92ZFRjlM6L0/y\n1b7mbUk+kuRxe7LtvU2S/ZO8OckN/d/CN5O8KklGvP+y/jlZOEv1zOrj3d8Z+nPj2VV1EPAkYCXw\n34c79AHQ+vE/AvhGVe2cYv29wMXAqVOsPwfYQffG8Xzg3UmO3sX2vgG8cKjtlL59PrwDOAN4OfBQ\n4NHAx4FnzdP258tHgFXAvwceDLwAOI1u/zVuVeW/WfwHXAf89sDy3wCf6G//E/BG4HPA3cBy4HBg\nA3AbsBV4ycB9zwTWAxcBPwKuAp4wsH4t8K1+3deA5w6se1G/nXcCdwBfB1YNrP8n4E8G+n52YF31\ntZ1GF7w7gDuBfwBeBXx0aJ//DnjHFMfj3/bb+gGwBXhO3/76/nHv7R/71F0c0+XdS/Xn2g7s7//o\ngbYLgLdM8RhnAuuAq4Gj+7aj++O2DjhzoO9L+ufitv65ObxvD3A28H3gh8BXgMdOdpwm2f5RwE+B\nY3axnwcDHwC2A9fTDRZ+aej5PLs/ltcCv96339jXdMrAY70feBfwyb6mzwEPB94O3N6/Hp443fM0\n8FjnAJfSvdY2Ab8yxT6sAu4Blg61H9vv//Ip/k7OBNb1t2+gew3e2f97KtO/nmf6eMuB/9M/1i3A\nRePOjvn61/pIc04lWUo32vniQPPEqOfBdH/YFwLb6ML/ROBNSZ450H813cjpocCHgI8n2a9f9y3g\nN+nC4vXAuiSPGLjvsX2fw4DXAR9L8tBR66+qc4EPAn9d3ZTPs+kC8vgkh/T7uBBYQxdWw/u/H90b\nxaeBhwEvAz6Y5DFV9TrgTXR/bAdV1Xmj1tV7NLCzqgZH6f9KF+S7cgH3jfZP6ZcHa34m8Gbgj4BH\ncN9zBPC7wNP7bR/c97l1iuM0bBWwrao+v4va/mf/uI8Cfquv88UD648FvgwcSvdauBB4Cl2A/THw\nzqGpxD+ie+M4DPgJcAXdwOEwusHE2/p9nvJ5GnisNXSvsYfQvSG+cYp9+B1gU1XdONhYVZvoXuer\ndrH/E57e//eQ/nheMbD/u/N6nuzx/opufx8CLKE79k0w9OfGx5P8APgs3WjiTQPr3l9VW6qb0ng4\n8DTgL6rqnqr6EvAefn4K4sqqWl9V99L9kR4AHAdQVR+pqpur6mdVdRHwTeCYgft+H3h7Vd3br7+G\nPZxKqKrvAP8MTMzBHw/cUlVXTtL9OOAgutH3jqq6HPgEcPKe1NA7iG60PegOujfTXVkHnNwH3Zp+\nedDzgfdW1VVV9RPgNcBTkyyjG80/GPhVIFV1dX88RnEoMGXfJAv6el5TVT+qquuA/0E3SJjw7ap6\nX1X9lO7T31LgrKr6SVV9mu6TxvKB/pdU1ZVVdQ9wCXBPVX1g4P5P7PuN8jxdUlWf71+3HwR+bYpd\nOWwX+/mdfv3ums3X871004uH9397U170sK8x9OfG71fVIVV1RFX9aVXdPbBucAR0OHBbVf1ooO16\nYPFk/avqZ9z3qYAkL0zypSQ/6N9kHsvP/1HdVP3n24HHPnyP9qxzPt3Ikv6/F0zR73Dgxr7uwRoW\nT9F/Ju4Efnmo7Zfpph+mVFU30I1U3wR8c3hESlfz9QP97wRuBRb3YfhOuqmO7yc5N8lwDVO5le6T\nw1QOA/Yb3Da/eKy+N3D77r6+4baDdtF/qr6jPE/fHbj946HtDLqFqffzEf363TWbr+dX003XfT7J\nliT/cQ/qul8x9Off4Iv2ZuChSQZHp48EbhpYXjpxoz/xuwS4OckRwN8Dp9Nd/XII8FW6F/KExUNX\nTDyy3+bu1jvh48DjkzwW+D26kd9kbgaWDp2wHt6/3fUNYGGSowbankA3Hz2dDwCvZJIpKbqaj5hY\nSHIg3Sj9JoCq+ruqejKwgm6a51V91+l+rvZ/A0uSrJxi/S3cN/qcMFvHajqz+Tx9Bji2n9r8/5Ic\nS/davrxvugt40ECXhw/cnupY7ur1PKPHq6rvVtVLqupw4KXAu5IsH+63LzL0x6gfZf5f4M1JDkjy\neLorVQanHJ6c5A/6ufM/p5ub/Re6E5lFd9KPJC+mG+kPehjw8iT7JTmJ7mTdxhmW+T26OebBuu+h\nmxP+EPD5fvQ8mU10o8JX9zU8A3g2982R71J/hdMBwAP65QOS7N/XcBfwMeCsJAcmeRrd+Y+pPnUM\nuohufv7iSdZ9GHhxkl/rt/Umujnq65I8Jcmx/dTQXXQnLCdGx79wnAZV1TfpTqx+OMkzkjyg3581\nSdb2Uy4XA29M8uD+Tf0V/OL001zYo+dpUFV9hu4N7qNJjk6yIMlxdPvx7v44AHwJWNNvbyXd+awJ\n2+mO6/Dx3NXreUaPl+SkJEv6xdvp/pYGP+nsswz98TsZWEY3YrkEeF3/hzPhfwH/ge6F+QLgD/o5\nza/RzfleQRc4j6O7umHQJrqrRm6hO/F2YlXdOsP6zgNW9FNIHx9oP7/f5pQhW1U76MLjhL6GdwEv\nrKqvj7jtI+imISZG73fTzeNO+FPggXRzvR8G/nNVTTvSr6q7q+ozQ9NuE+s+A/wl8FG6OehfoZtr\nh2766O/pnovr6aZs/qZfN9VxGvRy7pse+gHdScnn0p1Ehe4E6l10V+Z8lu5N9b3T7c+emoXnadgf\nAv8IfIpuGm4d3fF52UCfv6Q7trfTnSD+0EA9P6a/yq0/nsf1q3b1ep7p4z0F2JTkTrortM6o/rsp\n+7r8/BSZNJokj6S7bO7hVTV8QlWaVUleRHeJ8W+Mu5b7O0f6mrF+7vcVwIUGvnT/4teSNSP9ic3v\n0U1vHD/mciTNkNM7ktQQp3ckqSF73fTOYYcdVsuWLRt3GZJ0v3LllVfeUlWLpuu314X+smXL2Lx5\n87jLkKT7lSTXT9/L6R1JaoqhL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWrIXveN\n3D21bO2lc/r4171lj/6/4pI0Vo70Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x\n9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENf\nkhoyUugnOT7JNUm2Jlk7yfr9k1zUr9+UZFnfvl+S85N8JcnVSV4zu+VLkmZi2tBPsgA4BzgBWAGc\nnGTFULdTgdurajlwNvDWvv0kYP+qehzwZOClE28IkqT5N8pI/xhga1VdW1U7gAuB1UN9VgPn97fX\nA6uSBCjgwCQLgQcCO4AfzkrlkqQZGyX0FwM3Dixv69sm7VNVO4E7gEPp3gDuAr4D3AD8bVXdNryB\nJKcl2Zxk8/bt22e8E5Kk0cz1idxjgJ8ChwNHAq9M8qjhTlV1blWtrKqVixYtmuOSJKldo4T+TcDS\ngeUlfdukffqpnIOBW4HnAZ+qqnur6vvA54CVe1q0JGn3jBL6XwCOSnJkkgcAa4ANQ302AKf0t08E\nLq+qopvSeSZAkgOB44Cvz0bhkqSZmzb0+zn604HLgKuBi6tqS5Kzkjyn73YecGiSrcArgInLOs8B\nDkqyhe7N431V9eXZ3glJ0mgWjtKpqjYCG4faXjtw+x66yzOH73fnZO2SpPHwG7mS1BBDX5IaYuhL\nUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1\nxNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMM\nfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JashIoZ/k+CTXJNmaZO0k6/dP\nclG/flOSZQPrHp/kiiRbknwlyQGzV74kaSamDf0kC4BzgBOAFcDJSVYMdTsVuL2qlgNnA2/t77sQ\nWAf8p6o6GngGcO+sVS9JmpFRRvrHAFur6tqq2gFcCKwe6rMaOL+/vR5YlSTA7wJfrqp/BaiqW6vq\np7NTuiRppkYJ/cXAjQPL2/q2SftU1U7gDuBQ4NFAJbksyVVJXj3ZBpKclmRzks3bt2+f6T5IkkY0\n1ydyFwK/ATy//+9zk6wa7lRV51bVyqpauWjRojkuSZLaNUro3wQsHVhe0rdN2qefxz8YuJXuU8E/\nV9UtVfVjYCPwpD0tWpK0e0YJ/S8ARyU5MskDgDXAhqE+G4BT+tsnApdXVQGXAY9L8qD+zeC3gK/N\nTumSpJlaOF2HqtqZ5HS6AF8AvLeqtiQ5C9hcVRuA84ALkmwFbqN7Y6Cqbk/yNro3jgI2VtWlc7Qv\nkqRpTBv6AFW1kW5qZrDttQO37wFOmuK+6+gu25QkjZnfyJWkhhj6ktQQQ1+SGmLoS1JDDH1Jaoih\nL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS\n1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kN\nMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWrISKGf5Pgk1yTZmmTtJOv3T3JRv35TkmVD\n6x+Z5M4k/3V2ypYk7Y5pQz/JAuAc4ARgBXBykhVD3U4Fbq+q5cDZwFuH1r8N+OSelytJ2hOjjPSP\nAbZW1bVVtQO4EFg91Gc1cH5/ez2wKkkAkvw+8G1gy+yULEnaXaOE/mLgxoHlbX3bpH2qaidwB3Bo\nkoOAvwBev6sNJDktyeYkm7dv3z5q7ZKkGZrrE7lnAmdX1Z276lRV51bVyqpauWjRojkuSZLatXCE\nPjcBSweWl/Rtk/XZlmQhcDBwK3AscGKSvwYOAX6W5J6qeuceVy5JmrFRQv8LwFFJjqQL9zXA84b6\nbABOAa4ATgQur6oCfnOiQ5IzgTsNfEkan2lDv6p2JjkduAxYALy3qrYkOQvYXFUbgPOAC5JsBW6j\ne2OQJO1lRhnpU1UbgY1Dba8duH0PcNI0j3HmbtQnSZpFfiNXkhoy0khfo1m29tI5ffzr3vKsOX18\nSfs+R/qS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1J\naoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SG\nGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGjJS6Cc5Psk1SbYm\nWTvJ+v2TXNSv35RkWd/+O0muTPKV/r/PnN3yJUkzMW3oJ1kAnAOcAKwATk6yYqjbqcDtVbUcOBt4\na99+C/DsqnoccApwwWwVLkmauVFG+scAW6vq2qraAVwIrB7qsxo4v7+9HliVJFX1xaq6uW/fAjww\nyf6zUbgkaeZGCf3FwI0Dy9v6tkn7VNVO4A7g0KE+fwhcVVU/Gd5AktOSbE6yefv27aPWLkmaoXk5\nkZvkaLopn5dOtr6qzq2qlVW1ctGiRfNRkiQ1aZTQvwlYOrC8pG+btE+ShcDBwK398hLgEuCFVfWt\nPS1YkrT7Fo7Q5wvAUUmOpAv3NcDzhvpsoDtRewVwInB5VVWSQ4BLgbVV9bnZK1vDlq29dM63cd1b\nnjXn25A0t6Yd6fdz9KcDlwFXAxdX1ZYkZyV5Tt/tPODQJFuBVwATl3WeDiwHXpvkS/2/h836XkiS\nRjLKSJ+q2ghsHGp77cDte4CTJrnfG4A37GGN2suN81NGq9uWdpffyJWkhhj6ktQQQ1+SGjLSnL6k\nvYvnE7S7HOlLUkMMfUlqiNM7kmak1ctk95UpNUf6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGG\nviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhL\nUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSEjhX6S\n45Nck2RrkrWTrN8/yUX9+k1Jlg2se03ffk2Sfzd7pUuSZmra0E+yADgHOAFYAZycZMVQt1OB26tq\nOXA28Nb+viuANcDRwPHAu/rHkySNwSgj/WOArVV1bVXtAC4EVg/1WQ2c399eD6xKkr79wqr6SVV9\nG9jaP54kaQxSVbvukJwIHF9Vf9IvvwA4tqpOH+jz1b7Ptn75W8CxwJnAv1TVur79POCTVbV+aBun\nAaf1i48BrtnzXRvZYcAt87i9vYX73Rb3e993RFUtmq7TwvmoZDpVdS5w7ji2nWRzVa0cx7bHyf1u\ni/utCaNM79wELB1YXtK3TdonyULgYODWEe8rSZono4T+F4CjkhyZ5AF0J2Y3DPXZAJzS3z4RuLy6\neaMNwJr+6p4jgaOAz89O6ZKkmZp2eqeqdiY5HbgMWAC8t6q2JDkL2FxVG4DzgAuSbAVuo3tjoO93\nMfA1YCfwZ1X10znal901lmmlvYD73Rb3W8AIJ3IlSfsOv5ErSQ0x9CWpIc2G/nQ/LbGvSrI0yT8m\n+VqSLUnOGHdN8ynJgiRfTPKJcdcyX5IckmR9kq8nuTrJU8dd03xI8l/61/hXk3w4yQHjrmlv0GTo\nj/jTEvuqncArq2oFcBzwZw3tO8AZwNXjLmKevQP4VFX9KvAEGtj/JIuBlwMrq+qxdBehrBlvVXuH\nJkOf0X5aYp9UVd+pqqv62z+iC4DF461qfiRZAjwLeM+4a5kvSQ4Gnk53hR1VtaOqfjDequbNQuCB\n/XeHHgTcPOZ69gqthv5i4MaB5W00EnyD+l9DfSKwabyVzJu3A68GfjbuQubRkcB24H39tNZ7khw4\n7qLmWlXdBPwtcAPwHeCOqvr0eKvaO7Qa+s1LchDwUeDPq+qH465nriX5PeD7VXXluGuZZwuBJwHv\nrqonAncB+/w5rCQPofv0fiRwOHBgkj8eb1V7h1ZDv+mfh0iyH13gf7CqPjbueubJ04DnJLmObjrv\nmUnWjbekebEN2FZVE5/m1tO9Cezrfhv4dlVtr6p7gY8Bvz7mmvYKrYb+KD8tsU/qf/L6PODqqnrb\nuOuZL1X1mqpaUlXL6J7vy6tqnx/5VdV3gRuTPKZvWkX3Dfl93Q3AcUke1L/mV9HACexR7BW/sjnf\npvppiTGXNV+eBrwA+EqSL/Vt/62qNo6xJs2tlwEf7Ac41wIvHnM9c66qNiVZD1xFd8XaF/EnGQB/\nhkGSmtLq9I4kNcnQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ35f34Ocoio0vjMAAAAAElFTkSu\nQmCC\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAEICAYAAABbOlNNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAGb9JREFUeJzt3Xu4HHWd5/H3ZwgXBYUgDIQAgoIi\njgM6Z1gdL6sQlIsacEXRGcV9VNbdYVdHRwVxRnR1Jt4GXddxjKCgsAKikYwwIiDexssYNBICIhGj\n5AJERLwjyHf/qDradTiXJHVyOpD363n6OdVVv6rft7o7/en6VXUnVYUkSaP+aNgFSJI2LwaDJKnD\nYJAkdRgMkqQOg0GS1GEwSJI6DAbNuDQ+kuT2JP8x7HokdRkMm7kkK5P8OskvBm7/d9h19fQk4HBg\nz6o6ZOzCJHOSLE6yJkkl2WfM8m2TfDjJz5LcnOTVE3WU5CXtNk4fM39+O/+sPjuS5LQk56xHuxcm\nWdI+f2uT/FuSJ/Xpe3PTBv5rk9zQvmZ/lOQfk2y7AduoJPtNY03Tur0thcFw3/Csqtph4HbSdHeQ\nZNZ0b3MSDwVWVtUvJ1h+D/BZ4L9MsPw0YP92O08DXpfkiEn6+z7wvDH7eALwvQ0pemO1wfUe4B+A\n3YC9gX8G5s9E/zPo/wAnAi8GHgQcCRwGXDDMorQRqsrbZnwDVgLzJlj2EuArwLuA24EfAEcOLN8R\nOBNYC6wG3gpsNbDuvwOnA7eNLgPeDfy43dZJQAGzgOOAq8b0/2rgoglq2wNYDPwEWAG8vJ3/UuA3\nwO+AXwBvnmTfZ7X97zNm/hrg6QP3/zdw3hSP0WeBo9t5OwM3A+8Ezhpo+2xgOfBT4AvAowaWvb59\nDH8OXE/zhncE8FvgrnZfvjNO/zu2y46bZD+3pQmONe3tPcC27bKnAquA1wG3ts/lMcBRNMH2E+AN\nA9s6DfgEcE5b6zLgEcAp7fo3jXnsxn2eBrZ1AfDRdlvLgZEJ9mH/9jk9ZMz8vYA7gUPb+18AXjb2\n+Wmnv9Q+379sH7PnD+z/G2helyuBvxxYf0O3twvwmfY5/gnwZeCPhv3vfHO7ecRw3/efaN6odgHe\nAZyZJO2ys4C7gf2AxwJPB142Zt0baT7Fvg14Oc2nvIOBx9G8AY1aDOyb5FED815E86YxnvNo/kHv\nATwX+Ickh1bVmcArgK9Vc/Tzpg3Z2SSzgTnAdwZmfwd49BSrfpTmkyzA8cBFNG9Yo9t9BPBx4FXA\nrsAlwL8m2SbJI2lC8s+r6kHAM2iOeD5LcxRwfrsvB43T7xOA7YBFk9R2KvB4msf9IOAQ4I0Dy3dv\ntzEX+HvgQ8BfAX8GPBn4uyT7DrR/FvAxYDbwbeBSmtGBucBbgA8OtB33eRpY/uy2zU40r4GJhjEP\nA1ZVVeecUVXdBHydZuhwUlX1lHbyoPbxPH9g/3dp6z8BWNg+JxuzvdfQ7O+uNK/7N9CEhwYYDPcN\nn07y04HbyweW/bCqPlRVvwPOpnnT3C3JbjSfKl9VVb+sqltpjg6OH1h3TVW9r6rurqpfA88D3ltV\nq6rqdmDBaMOquhM4n+YNiSSPBvah+fTVkWQv4InA66vqN1W1FDiDP7wx97FD+/eOgXl30AxdTGYR\n8NQkO7Z1jA205wMXV9VlVXUXzVHYA4C/oPkkvC1wYJKtq2plVX1/Pet9CPDjqrp7kjZ/Cbylqm6t\nqnXAm2lCd9RdwNvaus6jeZN8b1X9vKqWA9fSBMqoL1fVpW2fn6B5E1wwsP4+SXZaz+fpK1V1Sfv6\n+tiYfgbtQnM0M5617fI+/q6q7qyqLwIX07xWN8ZdNP9GHlpVd1XVl6s9vNAfGAz3DcdU1U4Dtw8N\nLLt5dKKqftVO7kAz/r41sHY0UGg+Kf7xwLo3jelnjzHzxi4/G3hhe0TyIuCCNjDG2gP4SVX9fGDe\nD2k+8fX1i/bvgwfmPZhmqGNCbfBdTPNJ/CFV9e9jmuzR1jja/h6a/Z9bVStojiROA25Ncl6SPdaz\n3tuAXaY4h9Ppu50e3P5t7RszwK/bv7cMLP81fwjM8Zb9eJz1d2D9nqebB6Z/BWw3wb78mOYNdzxz\n2uUb6/bqno8a+/hsiHfSDJl9LsmNSU7uUdf9lsFw/3UTzVDJLgOB8uCqGhxyGftJaS2w58D9vQYX\nVtXXacbUnwy8kOYT5HjWADsnGfwUvzfNGH0v7ZHMWrqfXA+iGf+eykdphhLGu4poDU2YAs0VNjT7\nv7rt9/9V1ZPaNgW8fbSkKfr8Gs3zcMwkbTp90zxWa6bY7nSYzufp88BeSTpXmbVHJY8Hrmhn/RJ4\n4ECT3ddj27OTbD+mxtHHZ4O21x5lvaaqHkYzTPbqJIetRw1bFIPhfqqq1gKfA96d5MFJ/ijJw5P8\n50lWuwB4ZZK5SXaiOeE61kdpxpnvqqqvTND3TcBXgX9Msl2SP6U56TzlZZ2jkmxHM3wDsG17f7CG\nNyaZneQAmnMjZ63HZr9IM9b9vnGWXQAcneSwJFvTBMidwFeTPDLJoe1ll7+h+dR9T7veLTRDM+P+\nW6qqO2jOC7w/yTFJHphk6yRHJnlH2+zj7f7smmSXtv16P1Ybazqep4FtfQ/4F+DcJI9PslU73PhJ\n4PKqurxtuhR4Tvs47Nf2N+gW4GHjdPHm9nzPk4Fn0gyRbfD2kjwzyX5t8N9BM0x4D+owGO4b/jXd\n7zFMdiJz0IuBbWjGoG8HLmTiw31oTmp+Dria5qTlJTQnr3830OZjwJ8w9ZvHC2jOQayhGd9/08Cb\nw/r4NX8YNvoufxgCAXgTzSWoP6R5s39neyJ4UtW4oqp+Ms6y62nOn7yPZtjjWTSXCf+WJqAWtPNv\nphmOO6VddfQN6rYk35qg33fTXMH1RmAdzdHcScCn2yZvBZbQPO7LgG+182ZC3+dp0Ek05yjOoXnu\nPktz1dDgZcen0xx13kIzNHnumG2cBpzdDn+Onke4meb1u6Zt/4qq+u5Gbm9/4PK2vq8B/1xVV27k\n/t5vxfMumkiSI4F/qarBIZYH0Fz2+LiqumFoxWmLkOSpwDlVtedUbTV9PGLQ7yV5QJKjksxKMpfm\nk/nYo5P/DnzTUJDuv6YtGJIckeT6JCvGO9Of5mcMzm+XfyMDP3OQ5JR2/vVJnjFdNWmDheZSydtp\nhpKuoxnvbhYmK4FX0oy/S7qfmpahpCRb0XwL83CaL498E3hBVV070OZ/AH9aVa9IcjxwbFU9P8mB\nNCffDqG5BO1y4BEDl9dJkmbQdB0xHAKsqKob25N153Hv34GZT3NyCJqToIe1VwbMp/k5gzur6gc0\n1xjf64fVJEkzY7p+OG0u3S9DraL5uYVx21TV3UnuoPlW6Fyar8wPrnuvL0IlOZHmB7rYfvvt/+yA\nAw7YqEKXrb5j6kY9PWbujvZt30Pteyb6t++Z73uq/qdy1VVX/biqdp2q3Uz+omYvVbUQWAgwMjJS\nS5Ys2ajt7HPyxdNZ1riWLDjavu17qH3PRP/2PfN9T9X/VJL8cOpW0zeUtJrut2T35N7fnvx9m/Yr\n9TvS/FzA+qwrSZoh0xUM3wT2T7Jvkm1ofqht8Zg2i2l+GRGaX3H8fPvjVYuB49urlval+QKK/6uX\nJA3JtAwltecMTqL5ed+tgA9X1fIkbwGWVNVimv8X4GNJVtD8Dvrx7brLk1xA8+3cu4G/9ookSRqe\naTvHUFWX0PyEwuC8vx+Y/g3Nf/Yy3rpvo/n/ACRJQ+Y3nyVJHQaDJKnDYJAkdRgMkqQOg0GS1GEw\nSJI6DAZJUofBIEnqMBgkSR0GgySpw2CQJHUYDJKkDoNBktRhMEiSOgwGSVKHwSBJ6jAYJEkdBoMk\nqcNgkCR19AqGJDsnuSzJDe3f2eO0OTjJ15IsT3J1kucPLDsryQ+SLG1vB/epR5LUX98jhpOBK6pq\nf+CK9v5YvwJeXFWPBo4A3pNkp4Hlr62qg9vb0p71SJJ66hsM84Gz2+mzgWPGNqiq71XVDe30GuBW\nYNee/UqSNpG+wbBbVa1tp28GdpuscZJDgG2A7w/Mfls7xHR6km171iNJ6mnWVA2SXA7sPs6iUwfv\nVFUlqUm2Mwf4GHBCVd3Tzj6FJlC2ARYCrwfeMsH6JwInAuy9995TlS1J2khTBkNVzZtoWZJbksyp\nqrXtG/+tE7R7MHAxcGpVfX1g26NHG3cm+Qjwt5PUsZAmPBgZGZkwgCRJ/fQdSloMnNBOnwBcNLZB\nkm2ARcBHq+rCMcvmtH9Dc37imp71SJJ66hsMC4DDk9wAzGvvk2QkyRltm+cBTwFeMs5lqecmWQYs\nA3YB3tqzHklST1MOJU2mqm4DDhtn/hLgZe30OcA5E6x/aJ/+JUnTz28+S5I6DAZJUofBIEnqMBgk\nSR0GgySpw2CQJHUYDJKkDoNBktRhMEiSOgwGSVKHwSBJ6jAYJEkdBoMkqcNgkCR1GAySpA6DQZLU\nYTBIkjoMBklSh8EgSeowGCRJHb2DIcnOSS5LckP7d/YE7X6XZGl7Wzwwf98k30iyIsn5SbbpW5Mk\naeNNxxHDycAVVbU/cEV7fzy/rqqD29uzB+a/HTi9qvYDbgdeOg01SZI20nQEw3zg7Hb6bOCY9V0x\nSYBDgQs3Zn1J0vSbNQ3b2K2q1rbTNwO7TdBuuyRLgLuBBVX1aeAhwE+r6u62zSpg7ngrJzkROBFg\n7733noaypfu3lQuOHnYJW5z7y2O+XsGQ5HJg93EWnTp4p6oqSU2wmYdW1eokDwM+n2QZcMf6FlpV\nC4GFACMjIxP1IUnqab2CoarmTbQsyS1J5lTV2iRzgFsn2Mbq9u+NSb4APBb4JLBTklntUcOewOoN\n3AdJ0jSajnMMi4ET2ukTgIvGNkgyO8m27fQuwBOBa6uqgCuB5062viRp5kxHMCwADk9yAzCvvU+S\nkSRntG0eBSxJ8h2aIFhQVde2y14PvDrJCppzDmdOQ02SpI3U++RzVd0GHDbO/CXAy9rprwKPmWD9\nG4FD+tYhSZoefvNZktRhMEiSOgwGSVKHwSBJ6jAYJEkdBoMkqcNgkCR1GAySpA6DQZLUYTBIkjoM\nBklSh8EgSeowGCRJHQaDJKnDYJAkdRgMkqQOg0GS1GEwSJI6DAZJUkevYEiyc5LLktzQ/p09Tpun\nJVk6cPtNkmPaZWcl+cHAsoP71CNJ6q/vEcPJwBVVtT9wRXu/o6qurKqDq+pg4FDgV8DnBpq8dnR5\nVS3tWY8kqae+wTAfOLudPhs4Zor2zwX+rap+1bNfSdIm0jcYdquqte30zcBuU7Q/Hvj4mHlvS3J1\nktOTbDvRiklOTLIkyZJ169b1KFmSNJkpgyHJ5UmuGec2f7BdVRVQk2xnDvAY4NKB2acABwB/DuwM\nvH6i9atqYVWNVNXIrrvuOlXZkqSNNGuqBlU1b6JlSW5JMqeq1rZv/LdOsqnnAYuq6q6BbY8ebdyZ\n5CPA365n3ZKkTaTvUNJi4IR2+gTgoknavoAxw0htmJAkNOcnrulZjySpp77BsAA4PMkNwLz2PklG\nkpwx2ijJPsBewBfHrH9ukmXAMmAX4K0965Ek9TTlUNJkquo24LBx5i8BXjZwfyUwd5x2h/bpX5I0\n/fzmsySpw2CQJHUYDJKkjl7nGLRhVi44etglSDPC1/p9m0cMkqQOg0GS1GEwSJI6DAZJUofBIEnq\nMBgkSR0GgySpw2CQJHUYDJKkDoNBktRhMEiSOgwGSVKHP6K3hRjmj5ptqX1L91UeMUiSOgwGSVJH\n72BIclyS5UnuSTIySbsjklyfZEWSkwfm75vkG+3885Ns07cmSdLGm44jhmuA5wBfmqhBkq2A9wNH\nAgcCL0hyYLv47cDpVbUfcDvw0mmoSZK0kXoHQ1VdV1XXT9HsEGBFVd1YVb8FzgPmJwlwKHBh2+5s\n4Ji+NUmSNt5MnWOYC9w0cH9VO+8hwE+r6u4x8+8lyYlJliRZsm7duk1arCRtydbrctUklwO7j7Po\n1Kq6aHpLGl9VLQQWAoyMjNRM9ClJW6L1Coaqmtezn9XAXgP392zn3QbslGRWe9QwOl+SNCQzNZT0\nTWD/9gqkbYDjgcVVVcCVwHPbdicAM3IEIkka33RcrnpsklXAE4CLk1zazt8jySUA7dHAScClwHXA\nBVW1vN3E64FXJ1lBc87hzL41SZI2Xu+fxKiqRcCiceavAY4auH8JcMk47W6kuWpJkrQZ8JvPkqQO\ng0GS1GEwSJI6DAZJUofBIEnqMBgkSR0GgySpw2CQJHUYDJKkDoNBktRhMEiSOgwGSVKHwSBJ6jAY\nJEkdBoMkqcNgkCR1GAySpA6DQZLUYTBIkjp6BUOS45IsT3JPkpEJ2uyV5Mok17ZtXzmw7LQkq5Ms\nbW9HjbcNSdLMmdVz/WuA5wAfnKTN3cBrqupbSR4EXJXksqq6tl1+elW9q2cdkqRp0isYquo6gCST\ntVkLrG2nf57kOmAucO2EK0mShmZGzzEk2Qd4LPCNgdknJbk6yYeTzJ7JeiRJ9zZlMCS5PMk149zm\nb0hHSXYAPgm8qqp+1s7+APBw4GCao4p3T7L+iUmWJFmybt26DelakrQBphxKqqp5fTtJsjVNKJxb\nVZ8a2PYtA20+BHxmkjoWAgsBRkZGqm9NkqTxbfKhpDQnIM4ErquqfxqzbM7A3WNpTmZLkoao18nn\nJMcC7wN2BS5OsrSqnpFkD+CMqjoKeCLwImBZkqXtqm+oqkuAdyQ5GChgJfDf+tQjbU5WLjh62CVs\nkXzc++t7VdIiYNE489cAR7XTXwHGvWypql7Up39J0vTzm8+SpA6DQZLUYTBIkjoMBklSh8EgSeow\nGCRJHQaDJKnDYJAkdRgMkqQOg0GS1GEwSJI6DAZJUofBIEnqMBgkSR0GgySpw2CQJHUYDJKkDoNB\nktRhMEiSOgwGSVJHr2BIclyS5UnuSTIySbuVSZYlWZpkycD8nZNcluSG9u/sPvVIkvrre8RwDfAc\n4Evr0fZpVXVwVQ0GyMnAFVW1P3BFe1+SNES9gqGqrquq63tsYj5wdjt9NnBMn3okSf3N1DmGAj6X\n5KokJw7M362q1rbTNwO7TbSBJCcmWZJkybp16zZlrZK0RZs1VYMklwO7j7Po1Kq6aD37eVJVrU7y\nx8BlSb5bVZ3hp6qqJDXRBqpqIbAQYGRkZMJ2kqR+pgyGqprXt5OqWt3+vTXJIuAQmvMStySZU1Vr\nk8wBbu3blySpn00+lJRk+yQPGp0Gnk5z0hpgMXBCO30CsL5HIJKkTaTv5arHJlkFPAG4OMml7fw9\nklzSNtsN+EqS7wD/AVxcVZ9tly0ADk9yAzCvvS9JGqIph5ImU1WLgEXjzF8DHNVO3wgcNMH6twGH\n9alBkjS9/OazJKnDYJAkdRgMkqQOg0GS1GEwSJI6DAZJUofBIEnqMBgkSR0GgySpw2CQJHUYDJKk\nDoNBktRhMEiSOgwGSVKHwSBJ6jAYJEkdBoMkqcNgkCR1GAySpI5ewZDkuCTLk9yTZGSCNo9MsnTg\n9rMkr2qXnZZk9cCyo/rUI0nqb1bP9a8BngN8cKIGVXU9cDBAkq2A1cCigSanV9W7etYhSZomvYKh\nqq4DSLK+qxwGfL+qftinX0nSpjPT5xiOBz4+Zt5JSa5O8uEksydaMcmJSZYkWbJu3bpNW6UkbcGm\nDIYklye5Zpzb/A3pKMk2wLOBTwzM/gDwcJqhprXAuydav6oWVtVIVY3suuuuG9K1JGkDTDmUVFXz\npqmvI4FvVdUtA9v+/XSSDwGfmaa+JEkbaSaHkl7AmGGkJHMG7h5LczJbkjREfS9XPTbJKuAJwMVJ\nLm3n75HkkoF22wOHA58as4l3JFmW5GrgacDf9KlHktRf36uSFtG99HR0/hrgqIH7vwQeMk67F/Xp\nX5I0/fzmsySpw2CQJHUYDJKkjr4/iXGfs3LB0cMuQZI2ax4xSJI6DAZJUofBIEnqMBgkSR0GgySp\nw2CQJHUYDJKkDoNBktRhMEiSOgwGSVJHqmrYNWywJOuAH85gl7sAP57B/jYX7veWxf2+/3toVU35\nfyPfJ4NhpiVZUlUjw65jprnfWxb3W6McSpIkdRgMkqQOg2H9LBx2AUPifm9Z3G8BnmOQJI3hEYMk\nqcNgkCR1GAyTSHJEkuuTrEhy8rDrmQlJ9kpyZZJrkyxP8sph1zSTkmyV5NtJPjPsWmZKkp2SXJjk\nu0muS/KEYdc0U5L8Tfs6vybJx5NsN+yaNgcGwwSSbAW8HzgSOBB4QZIDh1vVjLgbeE1VHQg8Hvjr\nLWS/R70SuG7YRcyw9wKfraoDgIPYQvY/yVzgfwEjVfUnwFbA8cOtavNgMEzsEGBFVd1YVb8FzgPm\nD7mmTa6q1lbVt9rpn9O8ScwdblUzI8mewNHAGcOuZaYk2RF4CnAmQFX9tqp+OtyqZtQs4AFJZgEP\nBNYMuZ7NgsEwsbnATQP3V7GFvEGOSrIP8FjgG8OtZMa8B3gdcM+wC5lB+wLrgI+0Q2hnJNl+2EXN\nhKpaDbwL+BGwFrijqj433Ko2DwaDxpVkB+CTwKuq6mfDrmdTS/JM4NaqumrYtcywWcDjgA9U1WOB\nXwJbyvm02TSjAPsCewDbJ/mr4Va1eTAYJrYa2Gvg/p7tvPu9JFvThMK5VfWpYdczQ54IPDvJSpph\nw0OTnDPckmbEKmBVVY0eFV5IExRbgnnAD6pqXVXdBXwK+Ish17RZMBgm9k1g/yT7JtmG5qTU4iHX\ntMklCc1483VV9U/DrmemVNUpVbVnVe1D81x/vqru958eq+pm4KYkj2xnHQZcO8SSZtKPgMcneWD7\nuj+MLeTE+1RmDbuAzVVV3Z3kJOBSmqsVPlxVy4dc1kx4IvAiYFmSpe28N1TVJUOsSZvW/wTObT8A\n3Qj81yHXMyOq6htJLgS+RXM13rfx5zEAfxJDkjSGQ0mSpA6DQZLUYTBIkjoMBklSh8EgSeowGCRJ\nHQaDJKnj/wORAfQUewKmmAAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "stream",
          "text": [
            "Fraction of outputs displayed: 0.26\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4gj23-BnHEtS",
        "colab_type": "text"
      },
      "source": [
        "We see that, for a good choice of $\\gamma$ and $\\beta$, ground state is the most probable outcome.\n",
        "\n",
        "Try changing the values of $\\gamma$ and $\\beta$ away from the optimal ones. You'll see that this experiment no longer finds the ground state for us."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "iIUkW0w4FGfa",
        "colab_type": "text"
      },
      "source": [
        "### Exercise: Experiment with Different Numbers of Layers\n",
        "See if you can get a closer to the true ground state (i.e., a larger fraction of measurements yielding the minimal energy) by adding more layers to the circuit."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9JRyH-WDVfMp",
        "colab_type": "text"
      },
      "source": [
        "### Exercise: Try Ising Model on a different graph, or With Different Interaction Strengths\n",
        "Instead of a square lattice, you can try to formulate the Ising Model on any graph you like. This just changes which qubits you link in the $U(\\gamma, C)$ layer. Each edge of the graph could also come with a different interaction coefficient, so that instead of $\\exp(i\\pi \\gamma Z_iZ_j/2)$ for that edge you would have $\\exp(i\\pi \\gamma J_{ij}Z_iZ_j/2)$ for some matrix $J_{ij}$ of coefficients. Note that you have to change both the $U(\\gamma, C)$ layer and the definition of the energy function to make this work."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "--aUfkiaUb3S",
        "colab_type": "text"
      },
      "source": [
        "### Exercise: Repeat Using Sampling\n",
        "\n",
        "On real hardware we need to use sampling to estimate expectation values.\n",
        "\n",
        "Adjust your code so that sampling is used instead of wavefunction evaluation.\n",
        "\n",
        "How many samples do you need to take to get good results? Try different values.\n",
        "\n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "klGsAVyNtxl2",
        "colab_type": "text"
      },
      "source": [
        "### Exercise: Transverse field Ising Model\n",
        "The Ising Model with transverse field replaces the $\\sum h_i Z_i$ term with a $\\sum h_i X_i$ term. Can we use the QAOA here as well? What are the differences?\n",
        "This is no longer a classical problem: in general the ground state will now be a superposition of elements of the computational basis. Can you make a circuit that prepares a state close to the gound state?\n",
        "\n"
      ]
    }
  ]
}